OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
dmumps_driver.F File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine dmumps (id)
subroutine dmumps_set_infog (info, infog, comm, myid)
subroutine dmumps_print_icntl (id, lp)
subroutine dmumps_print_keep (id, lp)
subroutine dmumps_check_dense_rhs (idrhs, idinfo, idn, idnrhs, idlrhs)
subroutine dmumps_set_k221 (id)
subroutine dmumps_check_redrhs (id)

Function/Subroutine Documentation

◆ dmumps()

subroutine dmumps ( type (dmumps_struc) id)

Definition at line 19 of file dmumps_driver.F.

20 USE dmumps_ooc
21 USE mumps_memory_mod
22 USE dmumps_struc_def
23 USE dmumps_static_ptr_m ! For Schur pointer
24 USE dmumps_save_restore
25C
26!$ USE OMP_LIB
27C
28 IMPLICIT NONE
29C
30C =======
31C Purpose
32C =======
33C
34C TO SOLVE a SPARSE SYSTEM OF LINEAR EQUATIONS.
35C GIVEN AN UNSYMMETRIC, SYMMETRIC, OR SYMMETRIC POSITIVE DEFINITE
36C SPARSE MATRIX A AND AN N-VECTOR B, THIS SUBROUTINE SOLVES THE
37C SYSTEM A x = b or ATRANSPOSE x = b.
38C
39C List of main functionalities provided by the package:
40C ----------------------------------------------------
41C -Unsymmetric solver with partial pivoting (LU factorization)
42C -Symmetric positive definite solver (LDLT factorization)
43C -General symmetric solver with pivoting
44C -Either elemental or assembled matrix input
45C -Analysis/Factorization/Solve callable separately
46C -Deficient matrices (symmetric or unsymmetric)
47C -Rank revealing
48C -Null space basis computation
49C -Solution
50C -Return the Schur complement matrix while
51C also providing solution of interior problem
52C -Distributed input matrix and analysis phase
53C -Sequential or parallel MPI version (any number of processors)
54C -Error analysis and iterative refinement
55C -Out-of-Core factorization and solution
56C -Solution phase:
57C -Multiple Right-Hand-sides (RHS)
58C -Sparse RHS
59C -Distributed RHS
60C -Computation of selected entries of the inverse of
61C original matrix.
62C - Block Low-Rank (BLR) approximation based factorization
63C
64C Method
65C ------
66C The method used is a parallel direct method
67C based on a sparse multifrontal variant
68C of Gaussian elimination with partial numerical pivoting.
69C An initial ordering for the pivotal sequence
70C is chosen using the pattern of the matrix A + A^T and is
71C later modified for reasons of numerical stability. Thus this code
72C performs best on matrices whose pattern is symmetric, or nearly so.
73C For symmetric sparse matrices or for very unsymmetric and
74C very sparse matrices, other software might be more appropriate.
75C
76C
77C References :
78C -----------
79C
80C P. Amestoy, J.-Y. L'Excellent, G. Moreau, On exploiting sparsity of
81C multiple right-hand sides in sparse direct solvers,
82C SIAM Journal on Scientific Computing, volume 41, number 2,
83C pages A269-A291 (2019)
84C
85C G. Moreau, PhD Thesis, ENS-Lyon, University of Lyon,
86C On the solution phase of direct methods for sparse linear systems
87C with multiple sparse right-hand sides, December 10th, 2018
88C
89C P. Amestoy, A. Buttari, J.-Y. L'Excellent and T. Mary,
90C Performance and scalability of the block low-rank multifrontal
91C factorization on multicore architectures,
92C ACM Transactions on Mathematical Software (2018)
93C
94C T. Mary, PhD Thesis, University of Toulouse,
95C Block Low-Rank multifrontal solvers: complexity, performance, and
96C scalability, November 2017.
97C
98C S. de la Kethulle de Ryhove, P. Jaysaval and D.V. Shantsev,
99C P. R. Amestoy, J.-Y. L'Excellent and T. Mary,
100C Large-scale 3D EM modeling with a Block Low-Rank MUMPS solver,
101C Geophysical Journal International, volume 209, number 3,
102C pages 1558-1571 (2017) .
103C
104C P. Amestoy, A. Buttari, J.-Y. L'Excellent and T. Mary,
105C On the complexity of the Block Low-Rank multifrontal factorization,
106C SIAM Journal on Scientific Computing, volume 39,
107C number 4, pages A1710-A1740 (2017).
108C
109C P. Amestoy, R. Brossier, A. Buttari, J.-Y. L'Excellent, T. Mary,
110C L. Metivier, A. Miniussi, and S. Operto.
111C Fast 3D frequency-domain full waveform inversion with a parallel
112C Block Low-Rank multifrontal direct solver: application to OBC data
113C from the North Sea, Geophysics, 81(6):R363--R383, (2016).
114C
115C P. Amestoy, C. Ashcraft, O. Boiteau, A. Buttari, J.-Y. L'Excellent,
116C and C. Weisbecker.
117C Improving multifrontal methods by means of block low-rank representations.
118C SIAM Journal on Scientific Computing, 37(3):A1451--A1474 (2015).
119C
120C W. M. Sid-Lakhdar, PhD Thesis from Universite de Lyon prepared at ENS Lyon,
121C Scaling the solution of large sparse linear systems using multifrontal
122C methods on hybrid shared-distributed memory architectures (2014).
123C
124C P. Amestoy, J.-Y. L'Excellent, W. Sid-Lakhdar,
125C Characterizing asynchronous broadcast trees for multifrontal factorizations,
126C Workshop on Combinatorial Scientific Computing,
127C Lyon, France, July 21-23 (2014).
128C
129C P. Amestoy, J.-Y. L'Excellent, F.-H. Rouet, W. Sid-Lakhdar,
130C Modeling 1D distributed-memory dense kernels for an asynchronous
131C multifrontal sparse solver, High-Performance Computing for Computational
132C Science, VECPAR 2014, Eugene, Oregon, USA, June 30 - July 3 (2014).
133C
134C J.-Y. L'Excellent and W. M. Sid-Lakhdar,
135C Introduction of shared-memory parallelism in a distributed-memroy
136C multifrontal solver, Parallel Computing (40):3-4, pages 34-46 (2014).
137C
138C C. Weisbecker, PhD Thesis supported by EDF, INPT-IRIT,
139C Improving multifrontal solvers by means of algebraic block low-rank
140C representations (2013).
141C
142C E. Agullo, P. Amestoy, A. Buttari, A. Guermouche, G. Joslin, J.-Y.
143C L'Excellent, X. S. Li, A. Napov, F.-H. Rouet, M. Sid-Lakhdar, S. Wang, C.
144C Weisbecker, I. Yamazaki,
145C Recent Advances in Sparse Direct Solvers, 22nd Conference on Structural
146C Mechanics in Reactor Technology, San Francisco (2013).
147C
148C P. Amestoy, A. Buttari, G. Joslin, J.-Y. L'Excellent, W. Sid-Lakhdar, C.
149C Weisbecker, M. Forzan, C. Pozza, R. Perrin, V. Pellissier,
150C Shared memory parallelism and low-rank approximation techniques applied
151C applied to direct solvers in FEM simulation in <i>IEEE Transactions on
152C Magnetics</i>, IEEE, Special issue, Compumag 2013 (2013).
153C
154C L. Boucher, P. Amestoy, A, Buttari, F.-H. Rouet and M. Chauvin,
155C INTEGRAL/SPI data segmentation to retrieve sources intensity variations,
156C Astronomy & Astrophysics, Article 52, 20 pages,
157C http://dx.doi.org/10.1051/0004-6361/201219605 (2013).
158C
159C F.-H. Rouet, PhD thesis from INPT, Toulouse, France,
160C Memory and Performance issues in parallel multifrontal factorization and
161C triangular solutions with sparse right-hand sides (2014).
162C
163C J.-Y. L'Excellent, Habilitation thesis from ENS Lyon,
164C Multifrontal methods: Parallelism, Memory Usage and Numerical
165C Aspects (2012).
166C
167C P. Amestoy, I.S. Duff, J.-Y. L'Excellent, Y. Robert, F.H. Rouet
168C and B. Ucar, On computing inverse entries of a sparse matrix in
169C an out-of-core environment,
170C SIAM J. on Scientific Computing Vol. 34 N. 4, p. 1975-1999 (2012).
171C
172C Amestoy, Buttari, Duff, Guermouche, L'Excellent, and Ucar
173C The Multifrontal Method, Encyclopedia of Parallel Computing,
174C editor David Padua, Springer (2011).
175C
176C Amestoy, Buttari, Duff, Guermouche, L'Excellent, and Ucar
177C MUMPS, Encyclopedia of Parallel Computing,
178C editor David Padua, Springer (2011).
179C
180C Agullo, Guermouche and L'Excellent, Reducing the {I/O} Volume in
181C Sparse Out-of-core Multifrontal Methods}, SIAM SISC, Vol 31, Nb. 6,
182C 4774-4794 (2010).
183C
184C Amestoy, Duff, Guermouche, Slavova, Analysis of the Solution Phase of a
185C Parallel Multifrontal Approach, Parallel Computing, Vol. 36, 3--15 (2010).
186C
187C Tzvetomila Slavova, PhD from INPT prepared at CERFACS,
188C Parallel triangular solution in the out-of-core multifrontal approach
189C for solving large sparse linear systems, available as CERFACS
190C Report TH/PA/09/59 (2009).
191C
192C Agullo, Guermouche and L'Excellent, A Parallel Out-of-core Multifrontal
193C Method: Storage of Factors on Disk and Analysis of Models for an
194C Out-of-core Active Memory, Parallel Computing, Special Issue on Parallel
195C Matrix Algorithms, Vol. 34, Nb 6-8, 296--317 (2008).
196C
197C Emmanuel Agullo, PhD Thesis from LIP-Ecole Normale Superieure de Lyon,
198C On the Out-of-core Factorization of Large Sparse Matrices (Nov 2008).
199C
200C Amestoy, Duff, Ruiz, and Ucar, "A parallel
201C matrix scaling algorithm".
202C In proceedings of VECPAR'08-International Meeting-High
203C Performance Computing for Computational Science, (Jan 2008).
204C
205C Guermouche and L'Excellent, Constructing Memory-minimizing Schedules
206C for Multifrontal Methods, ACM TOMS, Vol. 32, Nb. 1, 17--32 (2006).
207C
208C Amestoy, Guermouche, L'Excellent, and Pralet,
209C Hybrid scheduling for the parallel solution
210C of linear systems. Vol 32 (2), pp 136-156 (2006).
211C
212C Stephane Pralet, PhD from INPT prepared at CERFACS,
213C Constrained orderings and scheduling for parallel sparse linear algebra,
214C available as CERFACS technical report, TH/PA/04/105, (Sept 2004).
215C
216C Abdou Guermouche, PhD Thesis from LIP-Ecole Normale Superieure de Lyon,
217C Etude et optimisation du comportement memoire dans les methodes paralleles
218C de factorisation de matrices creuses (2004).
219C
220C Guermouche, L'Excellent and Utard, Impact of Reordering on the Memory of a
221C Multifrontal Solver, Parallel Computing, Vol. 29, Nb. 9, 1191--1218 (2003).
222C
223C Amestoy, Duff, L'Excellent and Xiaoye S. Li, Impact of the Implementation
224C of MPI Point-to-Point Communications on the Performance of Two General
225C Sparse Solvers, Parallel Computing, Vol. 29, Nb 7, 833--847 (2003).
226C
227C Amestoy, Duff, L'Excellent and Xiaoye S. Li, Analysis and Comparison of
228C Two General Sparse Solvers for Distributed Memory Computers, ACM TOMS,
229C Vol. 27, Nb 4, 388--421 (2001).
230C
231C Amestoy, Duff, Koster and L'Excellent (2001),
232C A fully asynchronous multifrontal solver using distributed dynamic
233C scheduling, SIAM Journal of Matrix Analysis and Applications,
234C Vol 23, No 1, pp 15-41 (2001).
235C
236C Amestoy, Duff and L'Excellent (2000),
237C Multifrontal parallel distributed symmetric and unsymmetric solvers,
238C Comput. Methods in Appl. Mech. Eng., 184, 501-520 (2000)
239C
240C Amestoy, Duff and L'Excellent (1998),
241C Parallelisation de la factorisation LU de matrices
242C creuses non-symmetriques pour des architectures a memoire distribuee,
243C Calculateurs Paralleles Reseaux et systemes repartis,
244C Vol 10(5), 509-520 (1998).
245C
246C PARASOL Deliverable D2.1d (final report),
247C DMUMPS Version 3.1, A MUltifrontal Massively Parallel Solver,
248C PARASOL project, EU ESPRIT IV LTR project 20160, (June 1999).
249C
250C Jacko Koster, PhD from INPT prepared at CERFACS, On the parallel solution
251C and the reordering of unsymmetric sparse linear systems (1997).
252C
253C Vincent Espirat, Master's thesis from INPT(ENSEEIHT)-IRIT, Developpement
254C d'une approche multifrontale pour machines a memoire distribuee et
255C reseau heterogene de stations de travail (1996).
256C
257C Patrick Amestoy, PhD from INPT prepared at CERFACS, Factorization of large
258C sparse matrices based on a multifrontal approach in a multiprocessor
259C environment, Available as CERFACS report TH/PA/91/2 (1991).
260C
261C============================================
262C Argument lists and calling sequences
263C============================================
264C
265C There is only one entry:
266*
267* A Fortran 90 driver subroutine DMUMPS has been designed as a user
268* friendly interface to the multifrontal code.
269* This driver, in addition to providing the
270* normal functionality of a sparse solver, incorporates some
271* pre- and post-processing.
272* This driver enables the user to preprocess the matrix to obtain a
273* maximum
274* transversal so that the permuted matrix has a zero-free diagonal,
275* to perform prescaling
276* of the original matrix (a choice of scaling strategies is provided),
277* to use iterative refinement to improve the solution,
278* and finally to perform error analysis.
279*
280* The driver routine DMUMPS offers similar functionalities to other
281* sparse direct solvers, depending on the value of one of
282* its parameters (JOB). The main ones are:
283*
284* (i) JOB = -1
285C initializes an instance of the package. This must be
286C called before any other call to the package concerning that instance.
287C It sets default values for other
288C components of DMUMPS_STRUC, which may then be altered before
289C subsequent calls to DMUMPS.
290C Note that three components of the structure must always be set by the
291C user (on all processors) before a call with JOB=-1. These are
292C id%COMM,
293C id%SYM, and
294C id%PAR.
295C CNTL, ICNTL can then be modified (see documentation) by the user.
296C
297* A value of JOB = -1 cannot be combined with other values for JOB
298*
299* (ii) JOB = 1 accepts the pattern of matrix A and chooses pivots
300* from the diagonal using a selection criterion to
301* preserve sparsity. It uses the pattern of A + A^T
302* but ignores numerical values. It subsequently constructs subsidiary
303* information for the actual factorization by a call with JOB_=_2.
304* An option exists for the user to
305* input the pivot sequence, in which case only the necessary
306* information for a JOB = 2 entry will be generated. We call the JOB=1
307* entry, the analysis phase.
308C The following components of the structure define the centralized matrix
309C pattern and must be set by the user (on the host only)
310C before a call with JOB=1:
311C --- id%N, id%NZ (32-bit int) or id%NNZ (64-bit int),
312C id%IRN, and id%JCN
313C if the user wishes to input the structure of the
314C matrix in assembled format (ICNTL(5)=0, and ICNTL(18) $\neq$ 3),
315C --- id%ELTPTR, and id%ELTVAR
316C if the user wishes to input the matrix in elemental
317C format (ICNTL(5)=1).
318C A distributed matrix format is also available (see documentation)
319C
320* (iii) JOB = 2 factorizes a matrix A using the information
321* from a previous call with JOB = 1. The actual pivot sequence
322* used may differ slightly from that of this earlier call if A is not
323* diagonally dominant.
324*
325* (iv) JOB = 3 uses the factors generated by a JOB = 2 call to solve
326* a system of equations A X = B or A^T X =B, where X and B are matrices
327* that can be either dense or sparse.
328* The sparsity of B is exploited to limit the number of operations
329* performed during solution. When only part of the solution is
330* also needed (such as when computing selected entries of A^1) then
331* further reduction of the number of operations is performed.
332* This is particularly beneficial in the context of an
333* out-of-core factorization.
334*
335* (v) JOB = -2 frees all internal data allocated by the package.
336*
337* A call with JOB=3 must be preceded by a call with JOB=2,
338* which in turn must be preceded by a call with JOB=1, which
339* in turn must be preceded by a call with JOB=-1. Since the
340* information passed from one call to the next is not
341* corrupted by the second, several calls with JOB=2 for matrices
342* with the same sparsity pattern but different values may follow
343* a single call with JOB=1, and similarly several calls with JOB=3
344* can be used for different right-hand sides.
345* Values 4, 5, 6 for the parameter JOB can invoke combinations
346* of the three basic operations corresponding to JOB=1, 2 or 3.
347*
348C
349*********
350C --------------------------------------
351C Explicit interface needed for routines
352C using a target argument if they appear
353C in the same compilation unit.
354C --------------------------------------
355 INTERFACE
356 SUBROUTINE dmumps_check_dense_rhs
357 &(idrhs, idinfo, idn, idnrhs, idlrhs)
358 DOUBLE PRECISION, DIMENSION(:), POINTER :: idRHS
359 INTEGER, intent(in) :: idN, idNRHS, idLRHS
360 INTEGER, intent(inout) :: idINFO(:)
361 END SUBROUTINE dmumps_check_dense_rhs
362 SUBROUTINE dmumps_ana_driver( id )
363 USE dmumps_struc_def
364 TYPE (DMUMPS_STRUC), TARGET :: id
365 END SUBROUTINE dmumps_ana_driver
366 SUBROUTINE dmumps_fac_driver( id )
367 USE dmumps_struc_def
368 TYPE (DMUMPS_STRUC), TARGET :: id
369 END SUBROUTINE dmumps_fac_driver
370 SUBROUTINE dmumps_solve_driver( id )
371 USE dmumps_struc_def
372 TYPE (DMUMPS_STRUC), TARGET :: id
373 END SUBROUTINE dmumps_solve_driver
374 SUBROUTINE dmumps_print_icntl(id, LP)
375 USE dmumps_struc_def
376 TYPE (DMUMPS_STRUC), TARGET, INTENT(IN) :: id
377 INTEGER :: LP
378 END SUBROUTINE dmumps_print_icntl
379 END INTERFACE
380* MPI
381* ===
382 include 'mpif.h'
383 INTEGER MASTER
384 parameter( master = 0 )
385 INTEGER IERR
386*
387* ==========
388* Parameters
389* ==========
390 TYPE (DMUMPS_STRUC) :: id
391C
392C Main components of the structure are:
393C ------------------------------------
394C
395C (see documentation for a complete description)
396C
397C JOB is an INTEGER variable which must be set by the user to
398C characterize the factorization step. Possible values of JOB
399C are given below
400C
401C 1 Analysis: Ordering and symbolic factorization steps.
402C 2 Scaling and Numerical Factorization
403C 3 Solve and Error analysis
404C 4 Analysis followed by numerical factorization
405C 5 Numerical factorization followed by Solving step
406C 6 Analysis, Numerical factorization and Solve
407C
408C N is an INTEGER variable which must be set by the user to the
409C order n of the matrix A. It is not altered by the
410C subroutine.
411C
412C NZ / NNZ are INTEGER / INTEGER(8) variables which must be set by the user
413C to the number of entries being input, in case of centralized assembled
414C entry. It is not altered by the subroutine. Only used if
415C ICNTL(5).eq.0 and ICNTL(18) .ne. 3 (assembled matrix entry,
416C or, at least, centralized matrix graph during analysis).
417C
418C Restriction: NZ > 0 or NNZ > 0.
419C If NNZ is different from 0, NNZ is used. Otherwise, NZ is used.
420C
421C NELT is an INTEGER variable which must be set by the user to the
422C number of elements being input. It is not altered by the
423C subroutine. Only used if ICNTL(5).eq.1 (elemental matrix entry).
424C Restriction: NELT > 0.
425C
426C IRN and JCN are INTEGER arrays of length [N]NZ.
427C IRN(k) and JCN(k), k=1..[N]NZ must be set on entry to hold
428C the row and column indices respectively.
429C They are not altered by the subroutine except when ICNTL(6) = 1.
430C (in which case only the column indices are modified).
431C The arrays are only used if ICNTL(5).eq.0 (assembled entry)
432C or out-of-range.
433C
434C ELTPTR is an INTEGER array of length NELT+1.
435C ELTVAR is an INTEGER array of length ELTPTR(NELT+1)-1.
436C ELTPTR(I) points in ELTVAR to the first variable in the list of
437C variables that correspond to element I. ELTPTR(NELT+1) points
438C to the first unused location in ELTVAR.
439C The positions ELTVAR(I) .. ELTPTR(I+1)-1 contain the variables
440C for element I. No free space is allowed between variable lists.
441C ELTPTR/ELTVAR are not altered by the subroutine.
442C The arrays are only used if ICNTL(5).ne.0 (element entry).
443C
444C A is a DOUBLE PRECISION array of length [N]NZ.
445C The user must set A(k) to the value
446C of the entry in row IRN(k) and column JCN(k) of the matrix.
447C It is not altered by the subroutine.
448C (Note that the matrix can also be provided in a distributed
449C assembled input format)
450C
451C RHS is a DOUBLE PRECISION array of length N that is only accessed when
452C JOB = 3, 5, or 6. On entry, RHS(i)
453C must hold the i th component of the right-hand side of the
454C equations being solved.
455C On exit, RHS(i) will hold the i th component of the
456C solution vector. For other values of JOB, RHS is not accessed and
457C can be declared to have size one.
458C RHS should only be available on the host processor. If
459C it is associated on other processors, an error is raised.
460C (Note that the right-hand sides can also be provided in a
461C sparse format).
462C
463C COLSCA, ROWSCA are DOUBLE PRECISION
464C arrays of length N that are used to hold
465C the values used to scale the columns and the rows
466C of the original matrix, respectively.
467C These arrays need to be set by the user
468C only if ICNTL(8) is set to -1. If ICNTL(8)=0,
469C COLSCA and ROWSCA are not accessed and
470C so can be declared to have size one.
471C For any other values of ICNTL(8),
472C the scaling arrays are computed before
473C numerical factorization. The factors of the scaled matrix
474C diag(ROWSCA(i)) <A diag(COLSCA(i)) are computed.
475C
476C The workspace is automatically allocated by the package.
477C At the beginning of the numerical phase. If the user wants to increase
478C the allocated workspace (typically, numerical pivoting that leads to extra
479C storage, or previous call to MUMPS that failed because of
480C a lack of allocated memory),
481C we describe in the following how the user can modify the size
482C of the workspace:
483C 1/ The memory relaxation parameter
484C ICNTL(14) is designed to control the increase, with respect to the
485C estimations performed during analysis, in the size of the workspace
486C allocated during the numerical phase.
487C 2/ The user can also provide
488C a unique parameter, ICNTL(23), holding the maximum size of the total
489C workspace (in Megabytes) that the package is allowed to use internally.
490C In this case we try as much as possible to follow the indication given
491C by the relaxation parameter (ICNTL(14)).
492C
493C If ICNTL(23) is greater than 0
494C then MUMPS automatically computes the size of the internal working arrays
495C such that the storage for all MUMPS internal data is equal to ICNTL(23).
496C The relaxation ICNTL(14) is first applied to
497C the internal integer working array and communication buffer sizes;
498C the remaining available space is given to the real/complex
499C internal working arrays.
500C A lower bound of ICNTL(23) (if ICNTL(14) has not
501C been modified since the analysis) is given by INFOG(26).
502C
503C If ICNTL(23) is left to its default value 0
504C then each processor will allocate workspace based on
505C the estimates computed during the analysis (INFO(17)
506C if ICNTL(14) has not been modified since analysis,
507C or larger if ICNTL(14) was increased).
508C Note that these estimates are accurate in the sequential
509C version of {\tt MUMPS}, but that they can be inaccurate
510C in the parallel case. Therefore, in parallel, we recommend
511C to use ICNTL(23) and provide a value significantly larger
512C than INFOG(26).
513C --------------------------------------------------------------------------------
514C
515C CNTL is a DOUBLE PRECISION array of length 15
516C that contains control parameters and must be set by the user. Default
517C values for the components may be set by a call to DMUMPS(JOB=-1)
518C Details of the control parameters are given in DMUMPSID.
519C
520C ICNTL is an INTEGER array of length 60
521C that contains control parameters and must be set by the user. Default
522C values for the components may be set by a call to DMUMPS(JOB=-1)
523C Details of the control parameters are given in DMUMPSID.
524C
525C INFO is an INTEGER array of length 80 that need not be set by the
526C user. On return from DMUMPS, a value of zero for INFO(1)
527C indicates that the subroutine has performed successfully.
528C Details of the control parameters are given in DMUMPSID.
529C
530C RINFO is a DOUBLE PRECISION array of length 40 that need not be set by the
531C user. This array supplies information on the execution of DMUMPS.
532C Details of the control parameters are given in DMUMPSID.
533C
534C
535*
536*
537* ====================
538* .. Error Return ..
539* ====================
540*
541C MUMPS uses the following mechanism to process errors that
542C may occur during the parallel execution of the code.
543C If, during a call to MUMPS, an error occurs on a processor,
544C this processor informs all the other processors before they
545C return from the call.
546C In parts of the code where messages are sent asynchronously
547C (for example the factorization and solve phases),
548C the processor on which the error occurs sends a message
549C to the other processors with a specific error tag.
550C On the other hand, if the error occurs in a subroutine that
551C does not use asynchronous communication, the processor propagates
552C the error to the other processors.
553C On successful completion, a call to MUMPS will exit with the
554C parameter id%INFOG(1) set to zero.
555C A negative value for id%INFOG(1) indicates that an
556C error has been detected on one of the processors.
557C For example, if processor s returns with
558C INFO(1)= -8 and INFO(2)=1000, then processor s ran out of integer
559C workspace during the factorization and the size of the workspace
560C should be increased by 1000 at least.
561C The other processors are informed about this error and return with
562C INFO(1)=-1 (i.e., an error occurred on another processor) and
563C INFO(2)=s (i.e., the error occurred on processor s).
564C If several processors raised an error, those processors do not overwrite
565C INFO(1), i.e., only processors that did not produce an error will set
566C INFO(1) to -1 and INFO(2) to the rank of the processor having the most
567C negative error code.
568C
569C The behaviour is slightly different for the global information
570C parameters INFOG(1) and INFOG(2):
571C in the previous example, all processors would return with
572C INFOG(1)=-8 and INFOG(2)=1000.
573C
574C The possible error codes returned in INFO(1) (and INFOG(1))
575C are fully described in the documentation.
576C
577C A positive value of INFO(1) is associated with a warning message
578C which will be output on unit ICNTL(2) (see documentation).
579C
580C
581C .. Local variables ..
582C
583 INTEGER JOBMIN, JOBMAX, OLDJOB
584!$ INTEGER NOMP, NOMPMIN, NOMPMAX
585 INTEGER I, J, MP, LP, MPG, KEEP235SAVE, KEEP242SAVE,
586 & KEEP243SAVE, KEEP495SAVE, KEEP497SAVE
587 INTEGER(8) :: I8
588 LOGICAL LANA, LFACTO, LSOLVE, PROK, LPOK, FLAG, PROKG
589 LOGICAL NOERRORBEFOREPERM
590 LOGICAL UNS_PERM_DONE,I_AM_SLAVE
591C Saved communicator (pb of interference)
592 INTEGER COMM_SAVE
593C Local copy of JOB
594 INTEGER JOB
595 CHARACTER(LEN=20) :: FROM_C_INTERFACE_STRING
596 INTEGER, PARAMETER :: ICNTL18DIST_MIN = 1
597 INTEGER, PARAMETER :: ICNTL18DIST_MAX = 3
598 INTEGER, DIMENSION(:), ALLOCATABLE :: UNS_PERM_INV
599C TIMINGS
600 DOUBLE PRECISION TIMEG, TIMETOTAL
601 INTEGER(8) :: FILE_SIZE,STRUC_SIZE
602 INTEGER:: ICNTL16_LOC
603!$ INTEGER:: PREVIOUS_OMP_THREADS_NUM
604 IF (id%JOB .EQ. -200) THEN
605C -----------------
606C QUICK TERMINATION:
607C -----------------
608C In case of very serious error that cannot be recovered
609C (example MPI crash, signal to kill a process, etc.),
610C a user error handler may want to call MUMPS with
611C JOB=-200 to just suppress the existing OOC files (if any)
612C before terminating the application and killing all
613C processes. No checks are performed in this case,
614C data is not freed and MPI should not be called
615C (since MPI might no longer work)
616 CALL dmumps_ooc_clean_files(id,ierr)
617 RETURN
618 ENDIF
619 noerrorbeforeperm = .false.
620 uns_perm_done = .false.
621 lana = .false.
622 lfacto = .false.
623 lsolve = .false.
624 job = id%JOB
625C
626C Initialize error return codes to 0.
627 id%INFO(1) = 0
628 id%INFO(2) = 0
629C -----------------------------------
630C Check that MPI has been initialized
631C -----------------------------------
632 CALL mpi_initialized( flag, ierr )
633 IF ( .NOT. flag ) THEN
634 id%INFO(1) = -23
635 id%INFO(2) = 0
636 WRITE(6,990)
637 990 FORMAT(' Unrecoverable Error in DMUMPS initialization: ',
638 & ' MPI is not running.')
639 RETURN
640 END IF
641C ---------------------------
642C Duplicate user communicator
643C to avoid communications not
644C related to DMUMPS
645C ---------------------------
646 comm_save = id%COMM
647 CALL mpi_comm_dup( comm_save, id%COMM, ierr )
648C
649C Default setting for printing
650 lp = 6
651 mp = 0
652 mpg = 0
653 lpok = .true.
654 prok = .false.
655 prokg = .false.
656 icntl16_loc = 0
657C -------------------------
658C Check if value of JOB is
659C the same on all processes
660C -------------------------
661 CALL mpi_allreduce(job,jobmin,1,mpi_integer,mpi_max,
662 & id%COMM,ierr)
663 CALL mpi_allreduce(job,jobmax,1,mpi_integer,mpi_min,
664 & id%COMM,ierr)
665 IF ( jobmin .NE. jobmax ) THEN
666 id%INFO(1) = -3
667 id%INFO(2) = job
668 GOTO 499
669 END IF
670C
671C Check value of JOB and previous value of JOB
672C
673 IF ((job.LT.-3.OR.job.EQ.0.OR.job.GT.8)
674 & .AND. job.NE.9
675 & ) THEN
676C Out of range value
677 id%INFO(1) = -3
678 id%INFO(2) = job
679 GOTO 499
680 END IF
681 IF (job.NE.-1) THEN
682C Check the previous value of JOB
683C One should be able to test for old job value
684C Warning: non initialized value
685 oldjob = id%KEEP( 40 ) + 456789
686 IF (oldjob.NE.-1.AND.oldjob.NE.-2.AND.
687 & oldjob.NE.1.AND.oldjob.NE.2.AND.
688 & oldjob.NE.3) THEN
689 id%INFO(1) = -3
690 id%INFO(2) = job
691 GOTO 499
692 END IF
693 END IF
694 IF(job.NE.-1) THEN
695C Always allow JOB=-1
696 IF((job.GT.-2).AND.(id%KEEP(140).EQ.1)) then
697C If restore failed (keep(140)=1), the only allowed jobs (except -1)
698C are -2 or -3
699 id%INFO(1) = -3
700 id%INFO(2) = job
701 GOTO 499
702 END IF
703 ENDIF
704C ----------------------------------
705C Initialize, LANA, LFACTO, LSOLVE
706C LANA indicates if analysis must be performed
707C LFACTO indicates if factorization must be performed
708C LSOLVE indicates if solution must be performed
709C ----------------------------------
710 IF ((job.EQ.1).OR.(job.EQ.4).OR.
711 & (job.EQ.6)) lana = .true.
712 IF ((job.EQ.2).OR.(job.EQ.4).OR.
713 & (job.EQ.5).OR.(job.EQ.6)) lfacto = .true.
714 IF ((job.EQ.3).OR.(job.EQ.5).OR.
715 & (job.EQ.6)) lsolve = .true.
716 IF ( lana .OR. lfacto .OR. lsolve) THEN
717C Set some specific experimental KEEP entries,
718C Value defined on MASTER is the reference
719 CALL mpi_bcast( id%KEEP(370), 2, mpi_integer, master, id%COMM,
720 & ierr )
721 IF (lana) THEN
722 CALL mpi_bcast( id%KEEP(198), 1, mpi_integer, master,
723 & id%COMM, ierr )
724 IF (id%KEEP(370) .EQ. 1) THEN
725 id%KEEP(77)=0
726 id%KEEP(78)=0
727 id%KEEP(375)=1
728 ENDIF
729 ENDIF
730 IF (id%KEEP(371) .EQ. 1) THEN
731 IF (lana) THEN
732 IF (id%KEEP(50) .EQ. 0 .AND. id%NSLAVES .GE. 32) THEN
733 id%KEEP(376) = 1
734 ENDIF
735 ENDIF
736 IF (lfacto) THEN
737 id%KEEP(351) = 1
738 id%KEEP(206) = 2
739 ENDIF
740 IF (lsolve) THEN
741 id%KEEP(350) = 2
742 ENDIF
743 ENDIF
744 IF (lana) THEN
745 IF (id%KEEP(198).NE.0) THEN
746 id%KEEP(77)=0
747 id%KEEP(78)=0
748 id%KEEP(375)=1
749 IF ((id%KEEP(50).EQ.0) .AND. (id%NSLAVES.GT.1)) THEN
750 id%KEEP(376) = 1
751 ENDIF
752 IF (id%KEEP(198).EQ.2) THEN
753 id%KEEP(83) = 0
754 id%KEEP(91) = 0
755 ENDIF
756 ENDIF
757 ENDIF
758 ENDIF ! LANA .OR. LFACTO .OR. LSOLVE
759 IF (lana) THEN
760C Decode NZ / NNZ into KEEP8(28) now
761C since we may want to print it.
762C NZ and NNZ only accessed at analysis
763C phase
764 CALL mumps_get_nnz_internal( id%NNZ, id%NZ, id%KEEP8(28) )
765 ENDIF
766C Also decode NNZ_loc in the same way
767 IF (lana .OR. lfacto) THEN
768C NZ_loc/NNZ_loc may be accessed both analysis
769C and factorization:
770 CALL mumps_get_nnz_internal( id%NNZ_loc, id%NZ_loc,
771 & id%KEEP8(29))
772 ENDIF
773C
774 IF (job.EQ.-2.OR.job.EQ.1.OR.job.EQ.2.OR.job.EQ.3.OR.
775 & job.EQ.4.OR.job.EQ.5.OR.job.EQ.6
776 & .OR.job.EQ.7.OR.job.EQ.8.OR.job.EQ.-3
777 * .OR. job.EQ.9 ! LATER .OR.JOB.EQ.10
778 & ) THEN
779C Value of JOB is correct and differnet from -1
780C ICNTL should have been initialized and can be used
781 lp = id%ICNTL(1)
782 mp = id%ICNTL(2)
783 mpg = id%ICNTL(3)
784 lpok = ((lp.GT.0).AND.(id%ICNTL(4).GE.1))
785 prok = ((mp.GT.0).AND.(id%ICNTL(4).GE.2))
786 prokg = ( mpg .GT. 0 .and. id%MYID .eq. master )
787 prokg = (prokg.AND.(id%ICNTL(4).GE.2))
788 IF (id%KEEP(500).EQ.1) THEN
789 from_c_interface_string=" from C interface"
790 ELSE
791 from_c_interface_string=" "
792 ENDIF
793C -----------------------------
794C Set the number of threads if required.
795C -----------------------------
796 icntl16_loc = id%ICNTL(16)
797 CALL mpi_bcast( icntl16_loc, 1, mpi_integer, master, id%COMM,
798 & ierr )
799!$ IF (ICNTL16_LOC .GT. 0) THEN
800!$ PREVIOUS_OMP_THREADS_NUM = omp_get_max_threads()
801#if defined(WORKAROUNDINTELILP64OPENMPLIMITATION)
802!$ CALL omp_set_num_threads(int(ICNTL16_LOC,4))
803#else
804!$ CALL omp_set_num_threads(ICNTL16_LOC)
805#endif
806!$ ENDIF
807C -----------------------------
808C Check if number of threads is
809C the same on all processes
810C -----------------------------
811!$ NOMP = OMP_GET_MAX_THREADS()
812!$ CALL MPI_ALLREDUCE(NOMP,NOMPMIN,1,MPI_INTEGER,MPI_MIN,
813!$ & id%COMM,IERR)
814!$ CALL MPI_ALLREDUCE(NOMP,NOMPMAX,1,MPI_INTEGER,MPI_MAX,
815!$ & id%COMM,IERR)
816 id%KEEP(249) = 1
817!$ id%KEEP(249) = OMP_GET_MAX_THREADS()
818 IF (prokg) THEN
819C Print basic information on MUMPS call
820 IF (job .EQ. -2
821 & .OR. job .EQ. 8
822 & ) THEN
823C N, NELT, NNZ not meaningful
824C or not defined yet (JOB=8).
825 WRITE(mpg,'(/A,A,A,A,I4)')
826 & 'Entering DMUMPS ',
827 & trim(adjustl(id%VERSION_NUMBER)),
828 & trim(from_c_interface_string),
829 & ' with JOB =', job
830 ELSE IF (id%ICNTL(5) .NE. 1) THEN
831C Assembled format
832 IF (id%ICNTL(18) .EQ. 0
833 & ) THEN
834 WRITE(mpg,'(/A,A,A,A,I4,I12,I15)')
835 & 'Entering DMUMPS ',
836 & trim(adjustl(id%VERSION_NUMBER)),
837 & trim(from_c_interface_string),
838 & ' with JOB, N, NNZ =', job,id%N,id%KEEP8(28)
839 ELSE
840 WRITE(mpg,'(/A,A,A,A,I4,I12)')
841 & 'Entering DMUMPS ',
842 & trim(adjustl(id%VERSION_NUMBER)),
843 & trim(from_c_interface_string),
844 & ' with JOB, N =', job,id%N
845 ENDIF
846 ELSE
847C Elemental format
848 WRITE(mpg,'(/A,A,A,A,I4,I12,I15)')
849 & 'Entering DMUMPS ',
850 & trim(adjustl(id%VERSION_NUMBER)),
851 & trim(from_c_interface_string),
852 & ' driver with JOB, N, NELT =', job,id%N,id%NELT
853 ENDIF
854C MPI and OpenMP information
855!$ IF (.TRUE.) THEN
856!$ WRITE(MPG, '(A,I6,A,I6)') ' executing #MPI = ',
857!$ & id%NPROCS, ' and #OMP = ', NOMP
858!$ IF ( NOMPMIN .NE. NOMPMAX ) THEN
859!$ WRITE(MPG, '(A,I4,A,I4,A)')
860!$ & ' WARNING detected: different number of threads (max ',
861!$ & NOMPMAX, ', min ', NOMPMIN, ')'
862!$ END IF
863!$ ELSE
864 WRITE(mpg, '(A,I6,A)') ' executing #MPI = ',
865 & id%NPROCS, ', without OMP'
866!$ ENDIF
867 IF (job.GE.1 .AND. job.LE.6) THEN
868 IF ( id%KEEP(370).EQ.1.OR.id%KEEP(371).EQ.1) THEN
869 WRITE(mpg, 99996) id%KEEP(370), id%KEEP(371)
870 ENDIF
87199996 FORMAT(/'Advanced settings:'/
872 & ' KEEP(370) Static mapping =',i4/
873 & ' KEEP(371) Advanced optimizations =',i4)
874 IF (id%KEEP(401) .EQ.1) THEN
875 WRITE(mpg, 99997) id%KEEP(401)
876 ENDIF
87799997 FORMAT('L0 thread based multithreading setting:'/
878 & ' KEEP(401) (0=OFF, 1=ON) =',i4)
879 WRITE(mpg, '(A)')
880 ENDIF
881 ENDIF
882 END IF
883C
884C----------------------------------------------------------------
885C
886C JOB = -1 : START INITIALIZATION PHASE
887C (NEW INSTANCE)
888C
889C JOB = -2 : TERMINATE AN INSTANCE
890C----------------------------------------------------------------
891C
892 IF ( job .EQ. -1 ) THEN
893C
894C ------------------------------------------
895C Check that we have called (JOB=-2), ie
896C that the previous JOB is not 1 2 or 3,
897C before calling the initialization routine.
898C --------------------------------------------
899 id%INFO(1)=0
900 id%INFO(2)=0
901 oldjob = id%KEEP( 40 ) + 456789
902 IF ( oldjob .EQ. 1 .OR.
903 & oldjob .EQ. 2 .OR.
904 & oldjob .EQ. 3 ) THEN
905 IF ( id%N > 0 ) THEN
906 id%INFO(1)=-3
907 id%INFO(2)=job
908 ENDIF
909 ENDIF
910C Initialize id%MYID now because it is
911C required by MUMPS_PROPINFO. id%MYID
912C used to be initialized inside DMUMPS_INI_DRIVER,
913C leading to an uninitialized access here.
914 CALL mpi_comm_rank(id%COMM, id%MYID, ierr)
915 CALL mumps_propinfo( id%ICNTL(1),
916 & id%INFO(1),
917 & id%COMM, id%MYID )
918 IF ( id%INFO(1) .LT. 0 ) THEN
919C
920C If there was an error, then initialization
921C was already called and we can rely on the null
922C or non null value of the pointers related to OOC
923C stuff.
924C We use DMUMPS_CLEAN_OOC_DATA that should work even
925C on the master. Note that KEEP(201) was also
926C initialized in a previous call to Mumps.
927C
928C If DMUMPS_END_DRIVER or DMUMPS_FAC_DRIVER is called after
929C this error, then DMUMPS_CLEAN_OOC_DATA will be called
930C a second time, though.
931C
932 IF (id%KEEP(201).GT.0) THEN
933 CALL dmumps_clean_ooc_data(id, ierr)
934 ENDIF
935 GOTO 499
936 ENDIF
937C ----------------------------------------
938C Initialization DMUMPS_INI_DRIVER
939C ----------------------------------------
940C - Default values for ICNTL, KEEP,KEEP8, CNTL
941C - Attach emission buffer for buffered Send
942C - Nullify pointers in the structure
943C - Get rank and size of the communicator
944C ----------------------------------------
945 CALL dmumps_ini_driver( id )
946 IF ( id%INFO(1) .LT. 0 ) GOTO 499
947 GOTO 500
948 END IF
949 IF ( job .EQ. -2 ) THEN
950C -------------------------------------
951C Deallocation of the instance id
952C -------------------------------------
953 id%KEEP(40)= -2 - 456789
954 CALL dmumps_end_driver( id )
955 GOTO 500
956 END IF
957C
958C TIMINGS: for JOBS different from -1 and -2,
959C we measure TIMETOTAL:
960C
961 IF (id%MYID.EQ.master) THEN
962 id%DKEEP(70)=0.0d0
963 CALL mumps_secdeb(timetotal)
964 ENDIF
965C
966C----------------------------------------------------------------
967C
968C JOB = 7 : SAVE THE INSTANCE
969C
970C JOB = 8 : RESTORE THE INSTANCE
971C----------------------------------------------------------------
972C
973 IF ( job .EQ. 7 .OR. job .EQ. 8 ) THEN
974 IF( job.EQ.8 .AND. oldjob.NE.-1) THEN
975 id%INFO(1) = -3
976 id%INFO(2) = job
977 GOTO 499
978 END IF
979 IF (id%MYID.EQ.master) THEN
980C -----------------------------
981C Check incompatibility between
982C par (=0) and nprocs (=1)
983C -----------------------------
984 IF ( (id%KEEP(46).EQ.0).AND.(id%NPROCS.LE.1) )
985 & THEN
986 id%INFO(1) = -21
987 id%INFO(2) = id%NPROCS
988 ENDIF
989 ENDIF
990 CALL mumps_propinfo( id%ICNTL(1),
991 & id%INFO(1),
992 & id%COMM, id%MYID )
993 IF ( id%INFO(1) .LT. 0 ) GOTO 499
994 IF ( job .EQ. 7 ) THEN
995 IF (id%MYID.EQ.master) THEN
996 CALL mumps_secdeb(timeg)
997 ENDIF
998 CALL dmumps_save( id )
999 IF (id%MYID.EQ.master) THEN
1000 CALL mumps_secfin(timeg)
1001 IF (prokg) THEN
1002 WRITE( mpg,'(/A,F12.4)')
1003 & ' Elapsed time in save structure driver= ', timeg
1004 END IF
1005 ENDIF
1006 ELSE
1007 IF (id%MYID.EQ.master) THEN
1008 CALL mumps_secdeb(timeg)
1009 ENDIF
1010 CALL dmumps_restore( id )
1011 IF (id%MYID.EQ.master) THEN
1012 CALL mumps_secfin(timeg)
1013 IF (prokg) THEN
1014 WRITE( mpg,'(/A,F12.4)')
1015 & ' Elapsed time in restore structure driver= '
1016 & , timeg
1017 ENDIF
1018 END IF
1019 ENDIF
1020 IF ( id%INFO(1) .LT. 0 ) GOTO 499
1021 GOTO 500
1022 ENDIF
1023C
1024C----------------------------------------------------------------
1025C
1026C JOB = -3 : REMOVE SAVED INSTANCE
1027C
1028C----------------------------------------------------------------
1029C
1030 IF (job .EQ. -3) THEN
1031 CALL dmumps_remove_saved(id)
1032 IF ( id%INFO(1) .LT. 0 ) GOTO 499
1033 GOTO 500
1034 ENDIF
1035 IF (job.EQ.9) THEN
1036C Check that factorization was performed
1037 IF ( oldjob .LT. 2 ) THEN
1038 id%INFO(1)=-3
1039 id%INFO(2)=job
1040 ELSE
1041 CALL dmumps_sol_init_irhs_loc(id)
1042 ENDIF
1043 IF ( id%INFO(1) .LT. 0 ) GOTO 499
1044 GOTO 500
1045 ENDIF
1046C
1047C----------------------------------------------------------------
1048C
1049C MAIN DRIVER
1050C OTHER VALUES OF JOB : 1 to 6
1051C
1052C----------------------------------------------------------------
1053 CALL mumps_memory_set_data_sizes()
1054 IF (id%MYID.EQ.master) THEN
1055C -----------------------------
1056C Check incompatibility between
1057C par (=0) and nprocs (=1)
1058C -----------------------------
1059 IF ( (id%KEEP(46).EQ.0).AND.(id%NPROCS.LE.1) )
1060 & THEN
1061 id%INFO(1) = -21
1062 id%INFO(2) = id%NPROCS
1063 ENDIF
1064 END IF
1065C
1066C Propagate possible error to all nodes
1067 CALL mumps_propinfo( id%ICNTL(1),
1068 & id%INFO(1),
1069 & id%COMM, id%MYID )
1070 IF ( id%INFO(1) .LT. 0 ) GOTO 499
1071C
1072C Print ICNTL and KEEP
1073C
1074 IF (prok) CALL dmumps_print_icntl(id, mp)
1075C-----------------------------------------------------------------------
1076C
1077C CHECK SEQUENCE
1078C
1079C-----------------------------------------------------------------------
1080 IF ( lana ) THEN
1081 IF ( prokg .AND. oldjob .EQ. -1 ) THEN
1082C Print compilation options at first call to analysis
1083 CALL mumps_print_if_defined(mpg)
1084 ENDIF
1085C
1086C User wants to perform analysis. Previous value of
1087C JOB must be -1, 1, 2 or 3.
1088C
1089 IF ( oldjob .EQ. 0 .OR. oldjob .GT. 3 .OR. oldjob .LT. -1 ) THEN
1090 id%INFO(1) = -3
1091 id%INFO(2) = job
1092 GOTO 499
1093 END IF
1094 IF ( oldjob .GE. 2 ) THEN
1095C -----------------------------------------
1096C Previous step was factorization or solve.
1097C As analysis is now performed, deallocate
1098C at least some big arrays from facto.
1099C -----------------------------------------
1100 IF (associated(id%IS)) THEN
1101 DEALLOCATE (id%IS)
1102 NULLIFY (id%IS)
1103 END IF
1104 IF (associated(id%S)) THEN
1105 DEALLOCATE (id%S)
1106 NULLIFY (id%S)
1107 END IF
1108 END IF
1109 END IF
1110 IF ( lfacto ) THEN
1111C ------------------------------------
1112C User wants to perform factorization.
1113C Analysis must have been performed.
1114C ------------------------------------
1115 IF ( oldjob .LT. 1 .and. .NOT. lana ) THEN
1116 id%INFO(1) = -3
1117 id%INFO(2) = job
1118 GOTO 499
1119 END IF
1120 END IF
1121 IF ( lsolve ) THEN
1122C -------------------------------
1123C User wants to perform solve.
1124C Facto must have been performed.
1125C -------------------------------
1126 IF ( oldjob .LT. 2 .AND. .NOT. lfacto ) THEN
1127 id%INFO(1) = -3
1128 id%INFO(2) = job
1129 GOTO 499
1130 END IF
1131 END IF
1132C ------------------------------------------
1133C Permute JCN on entry to JOB if no analysis
1134C to be performed and IRN/JCN are needed.
1135C (facto: arrowheads + solve: iterative
1136C refinement and error analysis)
1137C ------------------------------------------
1138#if ! defined (LARGEMATRICES)
1139 noerrorbeforeperm =.true.
1140 uns_perm_done=.false.
1141 IF (id%MYID .eq. master .AND. id%KEEP(23) .NE. 0) THEN
1142 IF ( id%JOB .EQ. 2 .OR. id%JOB .EQ. 5 .OR.
1143 & (id%JOB .EQ. 3 .AND. (id%ICNTL(10) .NE.0 .OR.
1144 & id%ICNTL(11).NE. 0))) THEN
1145 uns_perm_done = .true.
1146 ALLOCATE(uns_perm_inv(id%N),stat=ierr)
1147 IF (ierr .GT. 0) THEN
1148C --------------------------------
1149C Exit with an error.
1150C We are not able to permute
1151C JCN correctly after a MAX-TRANS
1152C permutation resulting from a
1153C previous call to DMUMPS.
1154C --------------------------------
1155 id%INFO(1)=-13
1156 id%INFO(2)=id%N
1157 IF (lpok) WRITE(lp,99993)
1158 GOTO 510
1159 ENDIF
1160 DO i = 1, id%N
1161 uns_perm_inv(id%UNS_PERM(i))=i
1162 END DO
1163 DO i8 = 1_8, id%KEEP8(28)
1164 j = id%JCN(i8)
1165C -- skip out-of range (that are ignored in ANA_O)
1166 IF (j.LE.0.OR.j.GT.id%N) cycle
1167 id%JCN(i8)=uns_perm_inv(j)
1168 END DO
1169 DEALLOCATE(uns_perm_inv)
1170 END IF
1171 END IF
1172#endif
1173C
1174C Propagate possible error
1175 CALL mumps_propinfo( id%ICNTL(1),
1176 & id%INFO(1),
1177 & id%COMM, id%MYID )
1178 IF ( id%INFO( 1 ) .LT. 0 ) GO TO 499
1179*
1180*********
1181* MaxTrans-Analysis-Distri, Scale-Arrowhead-factorize, and
1182* Solve-IR-Error_Analysis (depending on the value of JOB)
1183*********
1184*
1185C
1186 IF ( lana ) THEN
1187C-----------------------------------------------------
1188C-
1189C- ANALYSIS : Max-Trans, Analysis, Distribution
1190C-
1191C-----------------------------------------------------
1192C
1193C Few checks + allocations
1194C
1195C IS : will be allocated on the slaves later
1196C PROCNODE : on the master only,
1197C because slave does not know N yet.
1198C Will be allocated in analysis for the slave.
1199C
1200C For assembled entry:
1201C IRN, JCN : check that they have been allocated by the
1202C user on the master, and if their size is adequate
1203C
1204C For element entry:
1205C ELTPTR, ELTVAR : check that they have been allocated by the
1206C user on the master, and if their size is adequate
1207C ----------------------------
1208C Reset KEEP(40) to -1 for the
1209C case where an error occurs
1210C ----------------------------
1211 id%KEEP(40)=-1 -456789
1212C
1213 IF (id%MYID.EQ.master) THEN
1214C Check N, [N]NZ, NELT
1215 IF ((id%N.LE.0).OR.((id%N+id%N+id%N)/3.NE.id%N)) THEN
1216 id%INFO(1) = -16
1217 id%INFO(2) = id%N
1218 GOTO 100
1219 END IF
1220 IF (id%ICNTL(5).NE.1) THEN
1221C Assembled input
1222 IF (id%ICNTL(18) .LT. 1 .OR. id%ICNTL(18) .GT. 3) THEN
1223C Centralized input
1224 IF (id%KEEP8(28) .LE. 0_8) THEN
1225 id%INFO(1) = -2
1226 CALL mumps_set_ierror(id%KEEP8(28), id%INFO(2))
1227 GOTO 100
1228 ENDIF
1229 ENDIF
1230 ELSE
1231C Element entry: check NELT on the master
1232 IF (id%NELT .LE. 0) THEN
1233 id%INFO(1) = -24
1234 id%INFO(2) = id%NELT
1235 GOTO 100
1236 ENDIF
1237 ENDIF
1238C -- initialize values of respectively
1239C icntl(6), (7) and (12) to not done/chosen
1240 id%INFOG(7) = -9999
1241 id%INFOG(23) = 0
1242 id%INFOG(24) = 1
1243C ---------------------------------------
1244C Element entry: allocate ELTPROC(1:NELT)
1245C ---------------------------------------
1246 IF ( id%ICNTL(5) .EQ. 1 ) THEN ! Elemental matrix
1247 IF ( associated( id%ELTPROC ) )
1248 & DEALLOCATE( id%ELTPROC )
1249 ALLOCATE( id%ELTPROC(id%NELT), stat=ierr )
1250 IF (ierr.gt.0) THEN
1251 id%INFO(1) = -7
1252 id%INFO(2) = id%NELT
1253 IF ( lpok ) WRITE(lp,'(A)')
1254 & 'Problem in allocating work array ELTPROC'
1255 GOTO 100
1256 END IF
1257 END IF
1258C ---------------------------------------------------
1259C Assembled centralized entry: check input parameters
1260C IRN/JCN
1261C Element entry: check input parameters ELTPTR/ELTVAR
1262C ---------------------------------------------------
1263 IF ( id%ICNTL(5) .NE. 1 ) THEN ! Assembled matrix
1264 id%KEEP8(30)=0_8
1265 IF ( id%ICNTL(18) .LT. icntl18dist_min
1266 & .OR. id%ICNTL(18) .GT. icntl18dist_max ) THEN
1267 IF ( .not. associated( id%IRN ) ) THEN
1268 id%INFO(1) = -22
1269 id%INFO(2) = 1
1270#if defined(MUMPS_F2003)
1271 ELSE IF ( size( id%IRN, kind=8 ) < id%KEEP8(28) ) THEN
1272#else
1273C size with kind=8 output not available before f2002. One can
1274C still check that if NZ can be stored in a 32-bit integer,
1275C the 32-bit size(id%IRN) is large enough
1276 ELSE IF ( id%KEEP8(28) .LE. int(huge(id%NZ),8) .AND.
1277 & size(id%IRN) < int(id%KEEP8(28)) ) THEN
1278#endif
1279 id%INFO(1) = -22
1280 id%INFO(2) = 1
1281 ELSE IF ( .not. associated( id%JCN ) ) THEN
1282 id%INFO(1) = -22
1283 id%INFO(2) = 2
1284#if defined(MUMPS_F2003)
1285 ELSE IF ( size( id%JCN, kind=8 ) < id%KEEP8(28) ) THEN
1286#else
1287C Same as for IRN above
1288 ELSE IF ( id%KEEP8(28) .LE. int(huge(id%NZ),8) .AND.
1289 & size(id%JCN) < int(id%KEEP8(28)) ) THEN
1290#endif
1291 id%INFO(1) = -22
1292 id%INFO(2) = 2
1293 END IF
1294 END IF
1295 IF ( id%INFO( 1 ) .eq. -22 ) THEN
1296 IF ( lpok ) WRITE(lp,'(A)')
1297 & 'Error in analysis: IRN/JCN badly allocated.'
1298 END IF
1299 ELSE
1300 IF ( .not. associated( id%ELTPTR ) ) THEN
1301 id%INFO(1) = -22
1302 id%INFO(2) = 1
1303 ELSE IF ( size( id%ELTPTR ) < id%NELT+1 ) THEN
1304 id%INFO(1) = -22
1305 id%INFO(2) = 1
1306 ELSE IF ( .not. associated( id%ELTVAR ) ) THEN
1307 id%INFO(1) = -22
1308 id%INFO(2) = 2
1309 ELSE
1310 id%LELTVAR = id%ELTPTR( id%NELT+1 ) - 1
1311 IF ( size( id%ELTVAR ) < id%LELTVAR ) THEN
1312 id%INFO(1) = -22
1313 id%INFO(2) = 2
1314 ELSE
1315C If no error, we compute KEEP8(30) (formerly NA_ELT),
1316C required for DMUMPS_MAX_MEM already in analysis, and
1317C then later during facto to check the size of A_ELT
1318 id%KEEP8(30) = 0_8
1319 IF ( id%KEEP(50) .EQ. 0 ) THEN
1320C Unsymmetric elements (but symmetric structure)
1321 DO i = 1,id%NELT
1322 j = id%ELTPTR(i+1) - id%ELTPTR(i)
1323 id%KEEP8(30) = id%KEEP8(30) + int(j,8) * int(j,8)
1324 ENDDO
1325 ELSE
1326C Symmetric elements
1327 DO i = 1,id%NELT
1328 j = id%ELTPTR(i+1) - id%ELTPTR(i)
1329 id%KEEP8(30) = id%KEEP8(30) +
1330 & (int(j,8) *int(j+1,8))/2_8
1331 ENDDO
1332 ENDIF
1333 ENDIF
1334 END IF
1335 IF ( id%INFO( 1 ) .eq. -22 ) THEN
1336 IF ( lpok ) WRITE(lp,'(A)')
1337 & 'Error in analysis: ELTPTR/ELTVAR badly allocated.'
1338 END IF
1339 ENDIF
1340 100 CONTINUE
1341 END IF
1342C
1343C Propagate possible error
1344 CALL mumps_propinfo( id%ICNTL(1),
1345 & id%INFO(1),
1346 & id%COMM, id%MYID )
1347 IF ( id%INFO( 1 ) .LT. 0 ) GO TO 499
1348C -----------------------------------------
1349C Call analysis procedure DMUMPS_ANA_DRIVER
1350C -----------------------------------------
1351 IF (id%MYID .eq. master) THEN
1352 id%DKEEP(71)=0.0d0
1353 CALL mumps_secdeb(timeg)
1354 END IF
1355C -------------------------------------------------
1356C Set scaling option for analysis in KEEP(52)
1357C (ICNTL(8) only defined on host at analysis phase)
1358C -------------------------------------------------
1359 IF (id%MYID.EQ.master) THEN
1360C{
1361 id%KEEP(52) = id%ICNTL(8)
1362C Out-of-range values => automatic choice
1363 IF ( id%KEEP(52) .GT. 8 .OR. id%KEEP(52).LT.-2)
1364 & id%KEEP(52) = 77
1365 IF ( id%KEEP(52) .EQ. 2 .OR. id%KEEP(52).EQ.5
1366 & .OR. id%KEEP(52) .EQ. 6 )
1367 & id%KEEP(52) = 77
1368 IF ((id%KEEP(52).EQ.77).AND.(id%KEEP(50).EQ.1)) THEN
1369 ! for SPD matrices default is no scaling
1370 id%KEEP(52) = 0
1371 ENDIF
1372 IF ( id%KEEP(52).EQ.77 .OR. id%KEEP(52).LE.-2) THEN
1373C -- suppress scaling computed during analysis
1374C -- if centralized matrix is not associated
1375 IF (.not.associated(id%A)) id%KEEP(52) = 0
1376 ENDIF
1377C deactivate analysis scaling if scaling given
1378 IF(id%KEEP(52) .EQ. -1) id%KEEP(52) = 0
1379C
1380C deactivate analysis scaling if
1381C permutation to zero-free diagonal not requested
1382 IF (id%ICNTL(6).EQ.0) id%KEEP(52) = 0
1383C deactivate analysis scaling for SPD matrices
1384 IF (id%KEEP(50).EQ.1) id%KEEP(52) = 0
1385C
1386 IF (id%KEEP(52).EQ.-2) THEN
1387C deallocate scalings in case of ordering allocated/computed
1388C during analysis. This is needed because in case of
1389C KEEP(52)=-2 then one cannot be sure that
1390C scaling will be effectivly computed during analysis
1391C Thus to test if scaling was effectively allocated/computed
1392C during analysis after DMUMPS_ANA_DRIVER one must
1393C be sure that scaling arrays are nullified.
1394 IF ( associated(id%COLSCA)) THEN
1395 DEALLOCATE( id%COLSCA )
1396 NULLIFY(id%COLSCA)
1397 ENDIF
1398 IF ( associated(id%ROWSCA)) THEN
1399 DEALLOCATE( id%ROWSCA )
1400 NULLIFY(id%ROWSCA)
1401 ENDIF
1402 ENDIF
1403C
1404C}
1405 ENDIF
1406C
1407C ANALYSIS PHASE:
1408 CALL dmumps_ana_driver( id )
1409C
1410C Check and save scaling option in INFOG(33)
1411 IF (id%MYID .eq. master) THEN
1412C{
1413 IF (id%KEEP(52).EQ.0) id%INFOG(33)=id%ICNTL(8)
1414 IF (id%KEEP(52).EQ.-2) THEN
1415C Scaling should have been computed during
1416 IF (.not.associated(id%COLSCA).OR.
1417 & .not.associated(id%ROWSCA)
1418 & ) THEN
1419C scaling was not computed reset KEEP(52)
1420C the user can then decide during factorization
1421C to activate scaling
1422 id%KEEP(52) =0
1423 id%INFOG(33)=0
1424 IF ( mpg .GT. 0 ) THEN
1425 WRITE(mpg,'(A)')
1426 & ' Warning; scaling was not computed during analysis'
1427 ENDIF
1428 IF ( associated(id%COLSCA)) THEN
1429 DEALLOCATE( id%COLSCA )
1430 NULLIFY(id%COLSCA)
1431 ENDIF
1432 IF ( associated(id%ROWSCA)) THEN
1433 DEALLOCATE( id%ROWSCA )
1434 NULLIFY(id%ROWSCA)
1435 ENDIF
1436 ENDIF
1437 ENDIF
1438 IF (id%KEEP(52) .NE. 0) THEN
1439 id%INFOG(33)=id%KEEP(52)
1440 ENDIF
1441C}
1442 ENDIF
1443C return value of ICNTL(12) effectively used
1444C that was saved on the master in KEEP(95)
1445 IF (id%MYID .eq. master) id%INFOG(24)=id%KEEP(95)
1446C TIMINGS:
1447 IF (id%MYID .eq. master) THEN
1448 CALL mumps_secfin(timeg)
1449 id%DKEEP(71) = timeg
1450 ENDIF
1451 IF (prokg) THEN
1452 WRITE( mpg,'(/A,F12.4)')
1453 & ' Elapsed time in analysis driver= ', timeg
1454 END IF
1455C -----------------------
1456C Return in case of error
1457C -----------------------
1458 IF ( id%INFO( 1 ) .LT. 0 ) GO TO 499
1459 id%KEEP(40) = 1 -456789
1460 END IF
1461C
1462C-------------------------------------------------------
1463C-
1464C
1465C BEGIN FACTORIZATION PHASE
1466C
1467C-
1468C-------------------------------------------------------
1469 IF ( lfacto ) THEN
1470 IF (id%MYID .eq. master) THEN
1471 id%DKEEP(91)=0.0d0
1472 CALL mumps_secdeb(timeg)
1473 END IF
1474C ----------------------
1475C Reset KEEP(40) to 1 in
1476C case of error in facto
1477C ----------------------
1478 id%KEEP(40) = 1 - 456789
1479C
1480C-------------------------------------------------------
1481C-
1482C- CHECKS, SCALING, ARROWHEAD + FACTORIZATION PHASE
1483C-
1484C-------------------------------------------------------
1485C
1486C Broadcast the value of KEEP(125) to decide if performing
1487C the scaling with the Schur complement feature.
1488 CALL mpi_bcast( id%KEEP(125), 1, mpi_integer, master, id%COMM,
1489 & ierr )
1490 IF ( id%MYID .EQ. master ) THEN
1491C -------------------------
1492C Check if Schur complement
1493C is allocated.
1494C -------------------------
1495 IF (id%KEEP(60).EQ.1) THEN
1496 IF ( associated( id%SCHUR_CINTERFACE)) THEN
1497C Called from C interface...
1498C Since id%SCHUR_CINTERFACE is of size 1,
1499C instruction below which causes bound check
1500C errors should be avoided. We cheat by first
1501C setting a static pointer with a routine with
1502C implicit interface, and then copying this pointer
1503C into id%SCHUR.
1504 CALL dmumps_set_tmp_ptr(id%SCHUR_CINTERFACE(1),
1505 & int(id%SIZE_SCHUR,8)*int(id%SIZE_SCHUR,8))
1506 CALL dmumps_get_tmp_ptr(id%SCHUR)
1507 NULLIFY(id%SCHUR_CINTERFACE)
1508 ENDIF
1509 IF ( .NOT. associated (id%SCHUR)) THEN
1510 IF (lp.GT.0)
1511 & write(lp,'(A)')
1512 & ' SCHUR not associated'
1513 id%INFO(1)=-22
1514 id%INFO(2)=9
1515 ELSE IF ( size(id%SCHUR) .LT.
1516 & id%SIZE_SCHUR * id%SIZE_SCHUR ) THEN
1517 IF (lp.GT.0)
1518 & write(lp,'(A)')
1519 & ' SCHUR allocated but too small'
1520 id%INFO(1)=-22
1521 id%INFO(2)=9
1522 END IF
1523 END IF
1524C ------------------------------------------------------------
1525C Assembled entry: check input parameterd IRN,JCN,A
1526C Element entry: check input parameters ELTPTR,ELTVAR,A_ELT
1527C ------------------------------------------------------------
1528 IF ( id%KEEP(54) .EQ. 0 ) THEN
1529 IF ( id%KEEP(55).eq.0 ) THEN
1530C Assembled entry
1531 IF ( .not. associated( id%IRN ) ) THEN
1532 id%INFO(1) = -22
1533 id%INFO(2) = 1
1534#if defined(MUMPS_F2003)
1535 ELSE IF ( size( id%IRN, kind=8 ) < id%KEEP8(28) ) THEN
1536#else
1537C size with kind=8 output not available. One can still
1538C check that if NZ can be stored in a 32-bit integer,
1539C the 32-bit size(id%IRN) (which we then assume not
1540C to overflow...) is large enough
1541 ELSE IF ( id%KEEP8(28) .LE. int(huge(id%NZ),8) .AND.
1542 & size(id%IRN) < int(id%KEEP8(28)) ) THEN
1543#endif
1544 id%INFO(1) = -22
1545 id%INFO(2) = 1
1546 ELSE IF ( .not. associated( id%JCN ) ) THEN
1547 id%INFO(1) = -22
1548 id%INFO(2) = 2
1549#if defined(MUMPS_F2003)
1550 ELSE IF ( size( id%JCN, kind=8 ) < id%KEEP8(28) ) THEN
1551#else
1552C Same as for IRN above
1553 ELSE IF ( id%KEEP8(28) .LE. int(huge(id%NZ),8) .AND.
1554 & size(id%JCN) < int(id%KEEP8(28)) ) THEN
1555#endif
1556 id%INFO(1) = -22
1557 id%INFO(2) = 2
1558 ELSEIF ( .not. associated( id%A ) ) THEN
1559 id%INFO( 1 ) = -22
1560 id%INFO( 2 ) = 4
1561#if defined(MUMPS_F2003)
1562 ELSE IF ( size( id%A, kind=8 ) < id%KEEP8(28) ) THEN
1563#else
1564C Same as for IRN/JCN above
1565 ELSE IF ( id%KEEP8(28) .LE. int(huge(id%NZ),8) .AND.
1566 & size( id%A ) < int(id%KEEP8(28)) ) THEN
1567#endif
1568 id%INFO( 1 ) = -22
1569 id%INFO( 2 ) = 4
1570 END IF
1571 ELSE
1572C Element entry
1573 IF ( .not. associated( id%ELTPTR ) ) THEN
1574 id%INFO(1) = -22
1575 id%INFO(2) = 1
1576 ELSE IF ( size( id%ELTPTR ) < id%NELT+1 ) THEN
1577 id%INFO(1) = -22
1578 id%INFO(2) = 1
1579 ELSE IF ( .not. associated( id%ELTVAR ) ) THEN
1580 id%INFO(1) = -22
1581 id%INFO(2) = 2
1582 ELSEIF ( size( id%ELTVAR ) < id%LELTVAR ) THEN
1583 id%INFO(1) = -22
1584 id%INFO(2) = 2
1585 ELSEIF ( .not. associated( id%A_ELT ) ) THEN
1586 id%INFO( 1 ) = -22
1587 id%INFO( 2 ) = 4
1588 ELSE
1589#if defined(MUMPS_F2003)
1590 IF ( size( id%A_ELT, kind=8 ) < id%KEEP8(30) ) THEN
1591#else
1592 IF ( id%KEEP8(30) < int(huge(id%NZ),8) .AND.
1593 & size( id%A_ELT ) < int(id%KEEP8(30)) ) THEN
1594#endif
1595 id%INFO( 1 ) = -22
1596 id%INFO( 2 ) = 4
1597 ENDIF
1598 END IF
1599 ENDIF
1600 ENDIF
1601C ----------------------
1602C Get the value of PERLU
1603C ----------------------
1604 CALL mumps_get_perlu(id%KEEP(12),id%ICNTL(14),
1605 & id%KEEP(50),id%KEEP(54),id%ICNTL(6),id%ICNTL(8))
1606C
1607C ----------------------
1608C Get null space options
1609C Note that nullspace is forbidden in case of Schur complement
1610C ----------------------
1611 CALL dmumps_get_ns_options_facto(id%N,id%KEEP(1),
1612 & id%ICNTL(1),mpg)
1613C ========================================
1614C Decode and set scaling options for facto
1615C ========================================
1616 IF (.NOT. ((id%KEEP(52).EQ.-2).AND.(id%ICNTL(8).EQ.77)) )
1617 & THEN
1618C if scaling was computed during analysis and automatic
1619C choice of scaling then we do not recompute scaling
1620 id%KEEP(52)=id%ICNTL(8)
1621 ENDIF
1622 IF ( id%KEEP(52) .GT. 8 .OR. id%KEEP(52).LT.-2)
1623 & id%KEEP(52) = 77
1624 IF ( id%KEEP(52) .EQ. 2 .OR. id%KEEP(52).EQ.5
1625 & .OR. id%KEEP(52) .EQ. 6 )
1626 & id%KEEP(52) = 77
1627 IF (id%KEEP(52).EQ.77) THEN
1628 IF (id%KEEP(50).EQ.1) THEN
1629 ! for SPD matrices the default is "no scaling"
1630 id%KEEP(52) = 0
1631 ELSE
1632 ! SYM .ne. 1 the default is cheap SIMSCA
1633 id%KEEP(52) = 7
1634 ENDIF
1635 ENDIF
1636 IF (id%KEEP(23) .NE. 0 .AND. id%ICNTL(8) .EQ. -1) THEN
1637 IF ( mpg .GT. 0 ) THEN
1638 WRITE(mpg,'(A)') ' ** WARNING : SCALING'
1639 WRITE(mpg,'(A)')
1640 & ' ** column permutation applied:'
1641 WRITE(mpg,'(A)')
1642 & ' ** column scaling has to be permuted'
1643 ENDIF
1644 ENDIF
1645C -----------------------------------
1646C If Schur has been asked for
1647C choose to disable or enable scaling
1648C ----------------------------------
1649 IF (id%KEEP(125).EQ.0) THEN
1650C ------------------------
1651C scaling is disabled
1652C ------------------------
1653 IF ( id%KEEP(60) .ne. 0 .and. id%KEEP(52) .ne. 0 ) THEN
1654 id%KEEP(52) = 0
1655 IF ( mpg .GT. 0 .AND. id%ICNTL(8) .NE. 0 ) THEN
1656 WRITE(mpg,'(A)') ' ** Warning: scaling not applied.'
1657 WRITE(mpg,'(A)') ' ** (disabled with Schur)'
1658 END IF
1659 END IF
1660 END IF
1661C -------------------------------
1662C If matrix is distributed on
1663C entry, only options 7 and 8
1664C of scaling are allowed.
1665C -------------------------------
1666 IF (id%KEEP(54) .NE. 0 .AND.
1667 & id%KEEP(52).NE.7 .AND. id%KEEP(52).NE.8 .AND.
1668 & id%KEEP(52) .NE. 0 ) THEN
1669 id%KEEP(52) = 0
1670 IF ( mpg .GT. 0 .and. id%ICNTL(8) .ne. 0 ) THEN
1671 WRITE(mpg,'(A)')
1672 & ' ** Warning: requested scaling option not available'
1673 WRITE(mpg,'(A)') ' ** for distributed matrix entry'
1674 END IF
1675 END IF
1676C ------------------------------------
1677C If matrix is symmetric, only scaling
1678C options -1 (given scaling), 1
1679C (diagonal scaling), 7 and 8 (SIMSCALING)
1680C are allowed.
1681C ------------------------------------
1682 IF ( id%KEEP(50) .NE. 0 ) THEN
1683 IF ( id%KEEP(52).ne. 1 .and.
1684 & id%KEEP(52).ne. -1 .and.
1685 & id%KEEP(52).ne. 0 .and.
1686 & id%KEEP(52).ne. 7 .and.
1687 & id%KEEP(52).ne. 8 .and.
1688 & id%KEEP(52).ne. -2 .and.
1689 & id%KEEP(52).ne. 77) THEN
1690 IF ( mpg .GT. 0 ) THEN
1691 WRITE(mpg,'(A)')
1692 & ' ** Warning: scaling option n.a. for symmetric matrix'
1693 END IF
1694 id%KEEP(52) = 0
1695 END IF
1696 END IF
1697C ----------------------------------
1698C If matrix is elemental on entry,
1699C automatic scaling is now forbidden
1700C ----------------------------------
1701 IF (id%KEEP(55) .NE. 0 .AND.
1702 & ( id%KEEP(52) .gt. 0 ) ) THEN
1703 id%KEEP(52) = 0
1704 IF ( mpg .GT. 0 ) THEN
1705 WRITE(mpg,'(A)') ' ** Warning: scaling not applied.'
1706 WRITE(mpg,'(A)')
1707 & ' ** (only user scaling av. for elt. entry)'
1708 END IF
1709 END IF
1710C --------------------------------------
1711C Check input parameters ROWSCA / COLSCA
1712C --------------------------------------
1713 IF ( id%KEEP(52) .eq. -1 ) THEN
1714 IF ( .not. associated( id%ROWSCA ) ) THEN
1715 id%INFO(1) = -22
1716 id%INFO(2) = 5
1717 ELSE IF ( size( id%ROWSCA ) < id%N ) THEN
1718 id%INFO(1) = -22
1719 id%INFO(2) = 5
1720 ELSE IF ( .not. associated( id%COLSCA ) ) THEN
1721 id%INFO(1) = -22
1722 id%INFO(2) = 6
1723 ELSE IF ( size( id%COLSCA ) < id%N ) THEN
1724 id%INFO(1) = -22
1725 id%INFO(2) = 6
1726 END IF
1727 END IF
1728C
1729C Allocate -- if required,
1730C ROWSCA and COLSCA on the master
1731C
1732C Allocation of scaling arrays.
1733C IF (KEEP(52)==-2 then scaling should have been allocated
1734C and computed during analysis
1735C
1736C If ICNTL(8) == -1, ROWSCA and COLSCA must have been associated and
1737C filled by the user. If ICNTL(8) is >0 and <= 8, the scaling is
1738C computed at the beginning of DMUMPS_FAC_DRIVER and is allocated now.
1739C
1740 IF (id%KEEP(52).GT.0 .AND.
1741 & id%KEEP(52) .LE.8) THEN
1742 IF ( associated(id%COLSCA))
1743 & DEALLOCATE( id%COLSCA )
1744 IF ( associated(id%ROWSCA))
1745 & DEALLOCATE( id%ROWSCA )
1746 ALLOCATE( id%COLSCA(id%N), stat=ierr)
1747 IF (ierr .GT.0) THEN
1748 id%INFO(1)=-13
1749 id%INFO(2)=id%N
1750 ENDIF
1751 ALLOCATE( id%ROWSCA(id%N), stat=ierr)
1752 IF (ierr .GT.0) THEN
1753 id%INFO(1)=-13
1754 id%INFO(2)=id%N
1755 ENDIF
1756 END IF
1757C
1758C Allocate scaling arrays of size 1 if
1759C they are not used to avoid problems
1760C when passing them in arguments
1761C
1762 IF (.NOT. associated(id%COLSCA)) THEN
1763 ALLOCATE( id%COLSCA(1), stat=ierr)
1764 END IF
1765 IF (ierr .GT.0) THEN
1766 id%INFO(1)=-13
1767 id%INFO(2)=1
1768 ENDIF
1769 IF (.NOT. associated(id%ROWSCA))
1770 & ALLOCATE( id%ROWSCA(1), stat=ierr)
1771 IF (ierr .GT.0) THEN
1772 id%INFO(1)=-13
1773 id%INFO(2)=1
1774 IF ( lpok ) WRITE(lp,'(A)')
1775 & 'problems in allocations before facto'
1776 GOTO 200
1777 END IF
1778.EQ. IF (id%KEEP(252) 1) THEN
1779 CALL DMUMPS_CHECK_DENSE_RHS
1780 & (id%RHS,id%INFO,id%N,id%NRHS,id%LRHS)
1781C Sets KEEP(221) and do some checks
1782 CALL DMUMPS_SET_K221(id)
1783 CALL DMUMPS_CHECK_REDRHS(id)
1784 ENDIF
1785 200 CONTINUE
1786.eq. END IF ! End of IF (MYID MASTER)
1787C KEEP(221) was set in DMUMPS_SET_K221 but not broadcast
1788 CALL MPI_BCAST( id%KEEP(221), 1, MPI_INTEGER, MASTER, id%COMM,
1789 & IERR )
1790C
1791C Check distributed matrices on all processors.
1792.ne..OR. I_AM_SLAVE = ( id%MYID MASTER
1793.eq..AND. & ( id%MYID MASTER
1794.eq. & id%KEEP(46) 1 ) )
1795.AND. IF (I_AM_SLAVE
1796.NE..AND..GT. & id%KEEP(54)0 id%KEEP8(29)0_8) THEN
1797.not. IF ( associated( id%IRN_loc ) ) THEN
1798 id%INFO(1) = -22
1799 id%INFO(2) = 16
1800#if defined(MUMPS_F2003)
1801 ELSE IF ( size( id%IRN_loc, KIND=8 ) < id%KEEP8(29) ) THEN
1802#else
1803C size with kind=8 output not available. One can still
1804C check that if NZ_loc can be stored in a 32-bit integer,
1805C the 32-bit size(id%IRN_loc) (which we then assume not
1806C to overflow...) is large enough
1807.LE..AND. ELSE IF ( id%KEEP8(29) int(huge(id%NZ_loc),8)
1808 & size(id%IRN_loc) < int(id%KEEP8(29)) ) THEN
1809#endif
1810 id%INFO(1) = -22
1811 id%INFO(2) = 16
1812.not. ELSE IF ( associated( id%JCN_loc ) ) THEN
1813 id%INFO(1) = -22
1814 id%INFO(2) = 16
1815#if defined(MUMPS_F2003)
1816 ELSE IF ( size( id%JCN_loc, KIND=8 ) < id%KEEP8(29) ) THEN
1817#else
1818C Same as for IRN_loc above
1819.LE..AND. ELSE IF ( id%KEEP8(29) int(huge(id%NZ_loc),8)
1820 & size(id%JCN_loc) < int(id%KEEP8(29)) ) THEN
1821#endif
1822 id%INFO(1) = -22
1823 id%INFO(2) = 16
1824.not. ELSEIF ( associated( id%A_loc ) ) THEN
1825 id%INFO( 1 ) = -22
1826 id%INFO( 2 ) = 16
1827#if defined(MUMPS_F2003)
1828 ELSE IF ( size( id%A_loc, KIND=8 ) < id%KEEP8(29) ) THEN
1829#else
1830C Same as for IRN_loc/JCN_loc above
1831.LE..AND. ELSE IF ( id%KEEP8(29) int(huge(id%NZ_loc),8)
1832 & size( id%A_loc ) < int(id%KEEP8(29)) ) THEN
1833#endif
1834 id%INFO( 1 ) = -22
1835 id%INFO( 2 ) = 16
1836 END IF
1837 ENDIF
1838C
1839C Check Schur complement on all processors.
1840C DMUMPS_PROPINFO will be called right after those checks.
1841C
1842.EQ..OR..EQ. IF (id%KEEP(60)2id%KEEP(60)3) THEN
1843 IF ( id%root%yes ) THEN
1844 IF ( associated( id%SCHUR_CINTERFACE )) THEN
1845C Called from C interface...
1846C The next instruction may cause
1847C bound check errors at runtime
1848C id%SCHUR=>id%SCHUR_CINTERFACE
1849C & (1:id%SCHUR_LLD*(id%root%SCHUR_NLOC-1)+
1850C & id%root%SCHUR_MLOC)
1851C Instead, we set a temporary
1852C pointer and then retrieve it
1853 CALL DMUMPS_SET_TMP_PTR(id%SCHUR_CINTERFACE(1),
1854 & int(id%SCHUR_LLD,8)*int(id%root%SCHUR_NLOC-1,8)+
1855 & int(id%root%SCHUR_MLOC,8))
1856 CALL DMUMPS_GET_TMP_PTR(id%SCHUR)
1857 NULLIFY(id%SCHUR_CINTERFACE)
1858 ENDIF
1859C Check that SCHUR_LLD is large enough
1860 IF (id%SCHUR_LLD < id%root%SCHUR_MLOC) THEN
1861.GT. IF (LP0) write(LP,*)
1862 & ' schur leading dimension schur_lld ',
1863 & id%SCHUR_LLD, 'too small with respect to',
1864 & id%root%SCHUR_MLOC
1865 id%INFO(1)=-30
1866 id%INFO(2)=id%SCHUR_LLD
1867.NOT. ELSE IF ( associated (id%SCHUR)) THEN
1868.GT. IF (LP0) write(LP,'(a)')
1869 & ' SCHUR not associated'
1870 id%INFO(1)=-22
1871 id%INFO(2)=9
1872 ELSE IF (size(id%SCHUR) <
1873 & id%SCHUR_LLD*(id%root%SCHUR_NLOC-1)+
1874 & id%root%SCHUR_MLOC) THEN
1875 IF (lp.GT.0) THEN
1876 write(lp,'(A)')
1877 & ' SCHUR allocated but too small'
1878 write(lp,*) id%MYID, ' : Size Schur=',
1879 & size(id%SCHUR),
1880 & ' SCHUR_LLD= ', id%SCHUR_LLD,
1881 & ' SCHUR_MLOC=', id%root%SCHUR_NLOC,
1882 & ' SCHUR_NLOC=', id%root%SCHUR_NLOC
1883 ENDIF
1884 id%INFO(1)=-22
1885 id%INFO(2)= 9
1886 ELSE
1887C We initialize the pointer that
1888C we will use within DMUMPS here.
1889 id%root%SCHUR_LLD=id%SCHUR_LLD
1890 IF (id%root%SCHUR_NLOC==0) THEN
1891 ALLOCATE(id%root%SCHUR_POINTER(1), stat=ierr)
1892 IF (ierr .GT.0) THEN
1893 id%INFO(1)=-13
1894 id%INFO(2)=1
1895 IF ( lpok ) THEN
1896 WRITE(lp,'(A)')
1897 & 'Problems in allocations before facto'
1898 ENDIF
1899 END IF
1900 ELSE
1901 id%root%SCHUR_POINTER=>id%SCHUR
1902 ENDIF
1903 ENDIF
1904 ENDIF
1905 ENDIF
1906C -------------------------
1907C Propagate possible errors
1908C -------------------------
1909 CALL mumps_propinfo( id%ICNTL(1),
1910 & id%INFO(1),
1911 & id%COMM, id%MYID )
1912 IF ( id%INFO(1) .LT. 0 ) GO TO 499
1913C -----------------------------------------------
1914C Call factorization procedure DMUMPS_FAC_DRIVER
1915C -----------------------------------------------
1916 CALL dmumps_fac_driver(id)
1917C Save scaling in INFOG(33)
1918 IF (id%MYID .eq. master) id%INFOG(33)=id%KEEP(52)
1919C
1920C In the case of Schur, free or not associated
1921C id%root%SCHUR_POINTER now rather than in end_driver.F
1922C (Case of repeated factorizations).
1923 IF (id%KEEP(60).EQ.2.OR.id%KEEP(60).EQ.3) THEN
1924 IF (id%root%yes) THEN
1925 IF (id%root%SCHUR_NLOC==0) THEN
1926 DEALLOCATE(id%root%SCHUR_POINTER)
1927 NULLIFY(id%root%SCHUR_POINTER)
1928 ELSE
1929 NULLIFY(id%root%SCHUR_POINTER)
1930 ENDIF
1931 ENDIF
1932 ENDIF
1933C root%RG2L_ROW and root%RG2L_COL
1934C are not used outside of the facto
1935 IF (associated(id%root%RG2L_ROW))THEN
1936 DEALLOCATE(id%root%RG2L_ROW)
1937 NULLIFY(id%root%RG2L_ROW)
1938 ENDIF
1939 IF (associated(id%root%RG2L_COL))THEN
1940 DEALLOCATE(id%root%RG2L_COL)
1941 NULLIFY(id%root%RG2L_COL)
1942 ENDIF
1943 IF (id%MYID .eq. master) THEN
1944 CALL mumps_secfin(timeg)
1945 id%DKEEP(91) = timeg
1946 ENDIF
1947 IF (prokg) THEN
1948 WRITE( mpg,'(/A,F12.4)')
1949 & ' Elapsed time in factorization driver= ', timeg
1950 END IF
1951C
1952C Check for errors after FACTO
1953C (it was propagated inside)
1954 IF(id%INFO(1).LT.0) THEN
1955C Free id%S if facto failed
1956 if (associated(id%S)) then
1957 DEALLOCATE(id%S)
1958 NULLIFY(id%S)
1959 endif
1960 GO TO 499
1961 ENDIF
1962C
1963C Update last successful step
1964C
1965 id%KEEP(40) = 2 - 456789
1966 END IF
1967C-------------------------------------------------------
1968C-
1969C
1970C BEGIN SOLVE PHASE
1971C
1972C-
1973C-------------------------------------------------------
1974 IF (lsolve) THEN
1975 IF (id%MYID .eq. master) THEN
1976 id%DKEEP(111)=0.0d0
1977 CALL mumps_secdeb(timeg)
1978 END IF
1979C ---------------------
1980C Reset KEEP(40) to 2.
1981C (last successful step
1982C was facto)
1983C ---------------------
1984 id%KEEP(40) = 2 -456789
1985C ------------------------------------------
1986C Call solution procedure DMUMPS_SOLVE_DRIVER
1987C ------------------------------------------
1988 IF (id%MYID .eq. master) THEN
1989 keep235save = id%KEEP(235)
1990 keep242save = id%KEEP(242)
1991 keep243save = id%KEEP(243)
1992 keep495save = id%KEEP(495)
1993 keep497save = id%KEEP(497)
1994 ! if no permutation of RHS asked then suppress request
1995 ! to interleave the RHS
1996 ! to interleave the RHS on ordering given then
1997 ! using option to set permutation to identity should be
1998 ! used (note though that
1999 ! they # with A-1/sparseRHS and Null Space)
2000 IF (id%KEEP(242).EQ.0) id%KEEP(243)=0
2001C --------------------------------------
2002C Check input parameters ROWSCA / COLSCA
2003C Only if KEEP(52).NE.0 because
2004C only 0 means that no colsca/rowsca are needed
2005C --------------------------------------
2006 IF ( id%KEEP(52) .ne. 0) THEN
2007 IF ( .not. associated( id%ROWSCA ) ) THEN
2008 id%INFO(1) = -22
2009 id%INFO(2) = 5
2010 ELSE IF ( size( id%ROWSCA ) < id%N ) THEN
2011 id%INFO(1) = -22
2012 id%INFO(2) = 5
2013 ELSE IF ( .not. associated( id%COLSCA ) ) THEN
2014 id%INFO(1) = -22
2015 id%INFO(2) = 6
2016 ELSE IF ( size( id%COLSCA ) < id%N ) THEN
2017 id%INFO(1) = -22
2018 id%INFO(2) = 6
2019 END IF
2020 ENDIF
2021 ENDIF
2022C -------------------------
2023C Propagate possible errors
2024C -------------------------
2025 CALL mumps_propinfo( id%ICNTL(1),
2026 & id%INFO(1),
2027 & id%COMM, id%MYID )
2028 IF ( id%INFO(1) .LT. 0 ) GO TO 499
2029 CALL dmumps_solve_driver(id)
2030 IF (id%MYID .eq. master) THEN
2031 CALL mumps_secfin(timeg)
2032 id%DKEEP(111) = timeg
2033 ENDIF
2034 IF (prokg) THEN
2035 WRITE( mpg,'(/A,F12.4)')
2036 & ' Elapsed time in solve driver= ', timeg
2037 END IF
2038 IF (id%MYID .eq. master) THEN
2039 id%KEEP(235) = keep235save
2040 id%KEEP(242) = keep242save
2041 id%KEEP(243) = keep243save
2042 id%KEEP(495) = keep495save
2043 id%KEEP(497) = keep497save
2044 ENDIF
2045 IF (id%INFO(1).LT.0) GOTO 499
2046C ---------------------------
2047C Update last successful step
2048C ---------------------------
2049 id%KEEP(40) = 3 -456789
2050 ENDIF
2051C
2052C What was actually done is saved in KEEP(40)
2053C
2054 IF (prok) CALL dmumps_print_icntl(id, mp)
2055 GOTO 500
2056*
2057*=================
2058* ERROR section
2059*=================
2060 499 CONTINUE
2061* Print error message if PROK
2062 IF (lpok) WRITE (lp,99995) id%INFO(1)
2063 IF (lpok) WRITE (lp,99994) id%INFO(2)
2064*
2065500 CONTINUE
2066#if ! defined(LARGEMATRICES)
2067C ---------------------------------
2068C Permute JCN on output to DMUMPS if
2069C KEEP(23) is different from 0.
2070C ---------------------------------
2071 IF (id%MYID .eq. master .AND. id%KEEP(23) .NE. 0
2072 & .AND. noerrorbeforeperm) THEN
2073C -------------------------------
2074C IF JOB=3 and PERM was not
2075C done (no iterative refinement/
2076C error analysis), then we do not
2077C permute JCN back.
2078C -------------------------------
2079 IF (id%JOB .NE. 3 .OR. uns_perm_done) THEN
2080 IF (.not.associated(id%UNS_PERM)) THEN
2081C I may happen
2082C (for ex in case of error -7 during analysis:
2083C UNS_PERM can be not associated,
2084C KEEP(23) was set to to automatic choice(=7) and
2085C an error of memory allocation occurs during analysis
2086C before having decided value of KEEP(23))
2087C UNS_PERM not associated and KEEP(23).NE.0
2088C Permuting JCN back does not make sense and KEEP(23)
2089C should be reset to zero
2090 id%KEEP(23) = 0
2091 ELSE
2092 DO i8 = 1_8, id%KEEP8(28)
2093 j=id%JCN(i8)
2094C -- skip out-of range (that are ignored in ANA_O)
2095 IF (j.LE.0.OR.j.GT.id%N) cycle
2096 id%JCN(i8)=id%UNS_PERM(j)
2097 END DO
2098 ENDIF
2099 END IF
2100 END IF
2101#endif
2102 510 CONTINUE
2103C ------------------------------------
2104C Set INFOG(1:2): same value on all
2105C processors + broadcast other entries
2106C ------------------------------------
2107 CALL dmumps_set_infog(id%INFO(1), id%INFOG(1), id%COMM, id%MYID)
2108C
2109C --------------------------------
2110C Broadcast RINFOG entries to make
2111C them available on all procs.
2112C --------------------------------
2113 CALL mpi_bcast( id%RINFOG(1), 40, mpi_double_precision, master,
2114 & id%COMM, ierr )
2115 IF (id%INFOG(1).GE.0 .AND. job.NE.-1
2116 & .AND. job.NE.-2 ) THEN
2117 IF (id%MYID .eq. master) THEN
2118 CALL mumps_secfin(timetotal)
2119 id%DKEEP(70) = timetotal
2120 ENDIF
2121 ENDIF
2122*=======================
2123* Compute space for save
2124*=======================
2125 IF (id%INFOG(1).GE.0) THEN
2126 CALL dmumps_compute_memory_save(id,file_size,struc_size)
2127 id%KEEP8(55)=file_size
2128 call mpi_allreduce(id%KEEP8(55),id%KEEP8(57),1,
2129 & mpi_integer8, mpi_sum,id%COMM,ierr)
2130 id%KEEP8(56)=struc_size
2131 call mpi_allreduce(id%KEEP8(56),id%KEEP8(58),1,
2132 & mpi_integer8, mpi_sum,id%COMM,ierr)
2133 id%RINFO(7)=dble(id%KEEP8(55))/1d6
2134 id%RINFO(8)=dble(id%KEEP8(56))/1d6
2135 id%RINFOG(17)=dble(id%KEEP8(57))/1d6
2136 id%RINFOG(18)=dble(id%KEEP8(58))/1d6
2137 ENDIF
2138!$ IF (ICNTL16_LOC .GT. 0) THEN
2139#if defined(WORKAROUNDINTELILP64OPENMPLIMITATION)
2140!$ CALL omp_set_num_threads(int(PREVIOUS_OMP_THREADS_NUM,4))
2141#else
2142!$ CALL omp_set_num_threads(PREVIOUS_OMP_THREADS_NUM)
2143#endif
2144!$ ICNTL16_LOC = 0
2145!$ ENDIF
2146*===============
2147* ERRORG section
2148*===============
2149 IF (id%MYID.EQ.master.and.mpg.GT.0.and.
2150 & id%INFOG(1).lt.0) THEN
2151 WRITE(mpg,'(A,I16)') ' On return from DMUMPS, INFOG(1)=',
2152 & id%INFOG(1)
2153 WRITE(mpg,'(A,I16)') ' On return from DMUMPS, INFOG(2)=',
2154 & id%INFOG(2)
2155 END IF
2156C -------------------------
2157C Restore user communicator
2158C -------------------------
2159 CALL mpi_comm_free( id%COMM, ierr )
2160 id%COMM = comm_save
2161 RETURN
2162*
216399995 FORMAT (' ** ERROR RETURN ** FROM DMUMPS INFO(1)=', i5)
216499994 FORMAT (' ** INFO(2)=', i16)
216599993 FORMAT (' ** Allocation error: could not permute JCN.')
mumps_propinfo
subroutine mumps_propinfo(icntl, info, comm, id)
Definition bcast_errors.F:15
dmumps_ana_driver
subroutine dmumps_ana_driver(id)
Definition dana_driver.F:16
dmumps_end_driver
subroutine dmumps_end_driver(id)
Definition dend_driver.F:15
dmumps_fac_driver
subroutine dmumps_fac_driver(id)
Definition dfac_driver.F:15
dmumps_ini_driver
subroutine dmumps_ini_driver(id)
Definition dini_driver.F:20
dmumps_set_infog
subroutine dmumps_set_infog(info, infog, comm, myid)
Definition dmumps_driver.F:2169
dmumps_print_icntl
subroutine dmumps_print_icntl(id, lp)
Definition dmumps_driver.F:2232
dmumps_check_dense_rhs
subroutine dmumps_check_dense_rhs(idrhs, idinfo, idn, idnrhs, idlrhs)
Definition dmumps_driver.F:2507
dmumps_get_ns_options_facto
subroutine dmumps_get_ns_options_facto(n, keep, icntl, mpg)
Definition drank_revealing.F:15
dmumps_sol_init_irhs_loc
subroutine dmumps_sol_init_irhs_loc(id)
Definition dsol_distrhs.F:505
dmumps_solve_driver
subroutine dmumps_solve_driver(id)
Definition dsol_driver.F:15
dmumps_set_tmp_ptr
subroutine dmumps_set_tmp_ptr(the_address, the_size8)
Definition dtools.F:1842
mpi_allreduce
subroutine mpi_allreduce(sendbuf, recvbuf, cnt, datatype, operation, comm, ierr)
Definition mpi.f:103
mpi_comm_dup
subroutine mpi_comm_dup(comm, comm2, ierr)
Definition mpi.f:230
mpi_bcast
subroutine mpi_bcast(buffer, cnt, datatype, root, comm, ierr)
Definition mpi.f:205
mpi_comm_free
subroutine mpi_comm_free(comm, ierr)
Definition mpi.f:238
mpi_comm_rank
subroutine mpi_comm_rank(comm, rank, ierr)
Definition mpi.f:254
mpi_initialized
subroutine mpi_initialized(flag, ierr)
Definition mpi.f:350
id
initmumps id
Definition multiplerhs_example.m:4
mumps_print_if_defined
subroutine mumps_print_if_defined(mpg)
Definition mumps_print_defined.F:15
dmumps_ooc
Definition dmumps_ooc.F:14
dmumps_ooc::dmumps_ooc_clean_files
subroutine dmumps_ooc_clean_files(id, ierr)
Definition dmumps_ooc.F:523
dmumps_ooc::dmumps_clean_ooc_data
subroutine dmumps_clean_ooc_data(id, ierr)
Definition dmumps_ooc.F:568
dmumps_save_restore
Definition dmumps_save_restore.F:14
dmumps_save_restore::dmumps_save
subroutine dmumps_save(id)
Definition dmumps_save_restore.F:300
dmumps_save_restore::dmumps_compute_memory_save
subroutine dmumps_compute_memory_save(id, total_file_size, total_struc_size)
Definition dmumps_save_restore.F:237
dmumps_save_restore::dmumps_remove_saved
subroutine dmumps_remove_saved(id)
Definition dmumps_save_restore.F:22
dmumps_save_restore::dmumps_restore
subroutine dmumps_restore(id)
Definition dmumps_save_restore.F:509
dmumps_static_ptr_m
Definition dstatic_ptr_m.F:14
dmumps_static_ptr_m::dmumps_get_tmp_ptr
subroutine, public dmumps_get_tmp_ptr(ptr)
Definition dstatic_ptr_m.F:24
dmumps_struc_def
Definition dmumps_struc_def.F:14
mumps_memory_mod
Definition mumps_memory_mod.F:14
mumps_memory_mod::mumps_memory_set_data_sizes
subroutine mumps_memory_set_data_sizes()
Definition mumps_memory_mod.F:26
mumps_secfin
subroutine mumps_secfin(t)
Definition tools_common.F:277
mumps_get_perlu
subroutine mumps_get_perlu(keep12, icntl14, keep50, keep54, icntl6, icntl8)
Definition tools_common.F:808
mumps_set_ierror
subroutine mumps_set_ierror(size8, ierror)
Definition tools_common.F:865
mumps_get_nnz_internal
subroutine mumps_get_nnz_internal(nnz, nz, nnz_i)
Definition tools_common.F:1248
mumps_secdeb
subroutine mumps_secdeb(t)
Definition tools_common.F:270
file_size
void file_size(int *filesize)
Definition write_routtines.c:371

◆ dmumps_check_dense_rhs()

subroutine dmumps_check_dense_rhs ( double precision, dimension(:), pointer idrhs,
integer, dimension(:), intent(inout) idinfo,
integer, intent(in) idn,
integer, intent(in) idnrhs,
integer, intent(in) idlrhs )

Definition at line 2505 of file dmumps_driver.F.

2507 IMPLICIT NONE
2508C
2509C Purpose:
2510C =======
2511C
2512C Check that the dense RHS is associated and of
2513C correct size. Called on master only, when dense
2514C RHS is supposed to be allocated. This can be used
2515C either at the beginning of the solve phase or
2516C at the beginning of the factorization phase
2517C if forward solve is done during factorization
2518C (see ICNTL(32)) ; idINFO(1), idINFO(2) may be
2519C modified.
2520C
2521C
2522C Arguments:
2523C =========
2524C
2525C id* : see corresponding components of the main
2526C MUMPS structure.
2527C
2528 DOUBLE PRECISION, DIMENSION(:), POINTER :: idRHS
2529 INTEGER, intent(in) :: idN, idNRHS, idLRHS
2530 INTEGER, intent(inout) :: idINFO(:)
2531 IF ( .not. associated( idrhs ) ) THEN
2532 idinfo( 1 ) = -22
2533 idinfo( 2 ) = 7
2534 ELSE IF (idnrhs.EQ.1) THEN
2535 IF ( size( idrhs ) < idn ) THEN
2536 idinfo( 1 ) = -22
2537 idinfo( 2 ) = 7
2538 ENDIF
2539 ELSE IF (idlrhs < idn)
2540 & THEN
2541 idinfo( 1 ) = -26
2542 idinfo( 2 ) = idlrhs
2543 ELSE IF
2544#if defined(MUMPS_F2003)
2545 & (size(idrhs,kind=8) <
2546 & int(idnrhs,8)*int(idlrhs,8)-idlrhs+idn)
2547#else
2548C size with kind=8 not available. One can still
2549C perform the check if minimal size small enough.
2550 & (int(idnrhs,8)*int(idlrhs,8)-idlrhs+idn
2551 & .LE. int(huge(idn),8)
2552 & .and.
2553 & size(idrhs) < int(int(idnrhs,8)*int(idlrhs,8)-idlrhs+idn))
2554#endif
2555 & THEN
2556 idinfo( 1 ) = -22
2557 idinfo( 2 ) = 7
2558 END IF
2559 RETURN

◆ dmumps_check_redrhs()

subroutine dmumps_check_redrhs ( type (dmumps_struc) id)

Definition at line 2584 of file dmumps_driver.F.

2585 USE dmumps_struc_def
2586 IMPLICIT NONE
2587C
2588C Purpose:
2589C =======
2590C
2591C * Decode API related to REDRHS and check REDRHS
2592C * Can be called at factorization or solve phase
2593C * Constraints:
2594C - Must be called after solve phase.
2595C - KEEP(60) must have been set (ok to check
2596C since KEEP(60) was set during analysis phase)
2597C * Remark that during solve phase, ICNTL(26)=1 is
2598C forbidden in case of fwd in facto.
2599C
2600 TYPE (DMUMPS_STRUC) :: id
2601 INTEGER MASTER
2602 parameter( master = 0 )
2603 IF (id%MYID .EQ. master) THEN
2604 IF ( id%KEEP(221) == 1 .or. id%KEEP(221) == 2 ) THEN
2605 IF (id%KEEP(221) == 2 .and. id%JOB == 2) THEN
2606 id%INFO(1)=-35
2607 id%INFO(2)=id%KEEP(221)
2608 GOTO 333
2609 ENDIF
2610 IF (id%KEEP(221) == 1 .and. id%KEEP(252) == 1
2611 & .and. id%JOB == 3) THEN
2612 id%INFO(1)=-35
2613 id%INFO(2)=id%KEEP(221)
2614 ENDIF
2615 IF ( id%KEEP(60).eq. 0 .or. id%SIZE_SCHUR.EQ.0 ) THEN
2616 id%INFO(1)=-33
2617 id%INFO(2)=id%KEEP(221)
2618 GOTO 333
2619 ENDIF
2620 IF ( .NOT. associated( id%REDRHS)) THEN
2621 id%INFO(1)=-22
2622 id%INFO(2)=15
2623 GOTO 333
2624 ELSE IF (id%NRHS.EQ.1) THEN
2625 IF (size(id%REDRHS) < id%SIZE_SCHUR ) THEN
2626 id%INFO(1)=-22
2627 id%INFO(2)=15
2628 GOTO 333
2629 ENDIF
2630 ELSE IF (id%LREDRHS < id%SIZE_SCHUR) THEN
2631 id%INFO(1)=-34
2632 id%INFO(2)=id%LREDRHS
2633 GOTO 333
2634 ELSE IF
2635 & (size(id%REDRHS)<
2636 & id%NRHS*id%LREDRHS-id%LREDRHS+id%SIZE_SCHUR)
2637 & THEN
2638 id%INFO(1)=-22
2639 id%INFO(2)=15
2640 GOTO 333
2641 ENDIF
2642 ENDIF
2643 ENDIF
2644 333 CONTINUE
2645C Error is not propagated. It should be propagated outside.
2646C The reason to propagate it outside is that there can be
2647C one call to PROPINFO instead of several ones.
2648 RETURN

◆ dmumps_print_icntl()

subroutine dmumps_print_icntl ( type (dmumps_struc), intent(in), target id,
integer lp )

Definition at line 2231 of file dmumps_driver.F.

2232 USE dmumps_struc_def
2233*
2234* Purpose:
2235* Print main control parameters CNTL and ICNTL
2236*
2237* ==========
2238* Parameters
2239* ==========
2240 TYPE (DMUMPS_STRUC), TARGET, INTENT(IN) :: id
2241 INTEGER :: LP
2242** Local Variables
2243 INTEGER, POINTER :: JOB
2244 INTEGER,DIMENSION(:),POINTER::ICNTL
2245 DOUBLE PRECISION, DIMENSION(:),POINTER::CNTL
2246 INTEGER MASTER
2247 parameter( master = 0 )
2248 IF (lp.LE.0) RETURN
2249 job=>id%JOB
2250 icntl=>id%ICNTL
2251 cntl=>id%CNTL
2252 IF (id%MYID.EQ.master) THEN
2253 SELECT CASE (job)
2254 CASE(1);
2255 WRITE (lp,980)
2256 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2257 IF (id%SYM.EQ.2) THEN
2258 WRITE (lp,991) icntl(5),icntl(6),icntl(7),icntl(12),
2259 & icntl(13),
2260 & icntl(15),
2261 & icntl(18),icntl(19),icntl(22),icntl(58)
2262 ELSE
2263 WRITE (lp,891) icntl(5),icntl(6),icntl(7),
2264 & icntl(13),
2265 & icntl(15),
2266 & icntl(18),icntl(19),icntl(22),icntl(58)
2267 ENDIF
2268 IF ((icntl(6).EQ.5).OR.(icntl(6).EQ.6).OR.
2269 & (icntl(12).NE.1) ) THEN
2270 WRITE (lp,992) icntl(8)
2271 ENDIF
2272 IF (id%ICNTL(19).NE.0)
2273 & WRITE(lp,998) id%SIZE_SCHUR
2274 WRITE (lp,993) icntl(14)
2275 CASE(2);
2276 WRITE (lp,980)
2277 WRITE (lp,981) cntl(1), cntl(3), cntl(4), cntl(5), cntl(7)
2278 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2279 WRITE (lp,992) icntl(8)
2280 WRITE (lp,993) icntl(14)
2281 WRITE (lp,923) icntl(24), icntl(31), icntl(32), icntl(33),
2282 & icntl(35), icntl(36)
2283 CASE(3);
2284 WRITE (lp,980)
2285 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2286 WRITE (lp,995)
2287 & icntl(9),icntl(10),icntl(11),icntl(20),icntl(21)
2288 CASE(4);
2289 WRITE (lp,980)
2290 WRITE (lp,981) cntl(1), cntl(3), cntl(4), cntl(5), cntl(7)
2291 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2292 WRITE (lp,992) icntl(8)
2293 IF (id%ICNTL(19).NE.0)
2294 & WRITE(lp,998) id%SIZE_SCHUR
2295 WRITE (lp,993) icntl(14)
2296 WRITE (lp,923) icntl(24), icntl(31), icntl(32), icntl(33),
2297 & icntl(35), icntl(36)
2298 CASE(5);
2299 WRITE (lp,980)
2300 WRITE (lp,981) cntl(1), cntl(3), cntl(4), cntl(5), cntl(7)
2301 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2302 IF (id%SYM.EQ.2) THEN
2303 WRITE (lp,991) icntl(5),icntl(6),icntl(7),icntl(12),
2304 & icntl(13),
2305 & icntl(15),
2306 & icntl(18),icntl(19),icntl(22),icntl(58)
2307 ELSE
2308 WRITE (lp,891) icntl(5),icntl(6),icntl(7),
2309 & icntl(13),
2310 & icntl(15),
2311 & icntl(18),icntl(19),icntl(22),icntl(58)
2312 ENDIF
2313 WRITE (lp,992) icntl(8)
2314 WRITE (lp,993) icntl(14)
2315 WRITE (lp,995)
2316 & icntl(9),icntl(10),icntl(11),icntl(20),icntl(21)
2317 WRITE (lp,923) icntl(24), icntl(31), icntl(32), icntl(33),
2318 & icntl(35), icntl(36)
2319 CASE(6);
2320 WRITE (lp,980)
2321 WRITE (lp,981) cntl(1), cntl(3), cntl(4), cntl(5), cntl(7)
2322 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2323 IF (id%SYM.EQ.2) THEN
2324 WRITE (lp,991) icntl(5),icntl(6),icntl(7),icntl(12),
2325 & icntl(13),
2326 & icntl(15),
2327 & icntl(18),icntl(19),icntl(22),icntl(58)
2328 ELSE
2329 WRITE (lp,891) icntl(5),icntl(6),icntl(7),
2330 & icntl(13),
2331 & icntl(15),
2332 & icntl(18),icntl(19),icntl(22),icntl(58)
2333 ENDIF
2334 IF (id%ICNTL(19).NE.0)
2335 & WRITE(lp,998) id%SIZE_SCHUR
2336 WRITE (lp,992) icntl(8)
2337 WRITE (lp,995)
2338 & icntl(9),icntl(10),icntl(11),icntl(20),icntl(21)
2339 WRITE (lp,993) icntl(14)
2340 WRITE (lp,923) icntl(24), icntl(31), icntl(32), icntl(33),
2341 & icntl(35), icntl(36)
2342 END SELECT
2343 ENDIF
2344 980 FORMAT (/'***********CONTROL PARAMETERS (ICNTL)**************'/)
2345 981 FORMAT (
2346 & ' CNTL(1) Threshold for numerical pivoting =',d16.4/
2347 & ' CNTL(3) Null pivot detection threshold =',d16.4/
2348 & ' CNTL(4) Threshold for static pivoting =',d16.4/
2349 & ' CNTL(5) Fixation for null pivots =',d16.4/
2350 & ' CNTL(7) Dropping threshold for BLR compression =',d16.4)
2351 990 FORMAT (
2352 & 'ICNTL(1) Output stream for error messages =',i10/
2353 & 'ICNTL(2) Output stream for diagnostic messages =',i10/
2354 & 'ICNTL(3) Output stream for global information =',i10/
2355 & 'ICNTL(4) Level of printing =',i10)
2356 991 FORMAT (
2357 & 'ICNTL(5) Matrix format ( keep(55) ) =',i10/
2358 & 'ICNTL(6) Maximum transversal ( keep(23) ) =',i10/
2359 & 'ICNTL(7) Ordering =',i10/
2360 & 'ICNTL(12) LDLT ordering strat ( keep(95) ) =',i10/
2361 & 'ICNTL(13) Parallel root (0=on, 1=off) =',i10/
2362 & 'ICNTL(15) Analysis by block =',i10/
2363 & 'ICNTL(18) Distributed matrix ( keep(54) ) =',i10/
2364 & 'ICNTL(19) Schur option ( keep(60) 0=off,else=on ) =',i10/
2365 & 'icntl(22) out-of-core option(0=off, >0=on) =',I10/
2366 & 'icntl(58) symbolic factorization option =',I10)
2367 891 FORMAT (
2368 & 'icntl(5) matrix format ( keep(55) ) =',I10/
2369 & 'icntl(6) maximum transversal( keep(23) ) =',I10/
2370 & 'icntl(7) ordering =',I10/
2371 & 'icntl(13) parallel root(0=on, 1=off) =',I10/
2372 & 'icntl(15) analysis by block =',I10/
2373 & 'icntl(18) distributed matrix( keep(54) ) =',I10/
2374 & 'icntl(19) schur option( keep(60) 0=off,else=on ) =',I10/
2375 & 'icntl(22) out-of-core option(0=off, >0=on) =',I10/
2376 & 'icntl(58) symbolic factorization option =',I10)
2377 992 FORMAT (
2378 & 'icntl(8) scaling strategy =',I10)
2379 923 FORMAT (
2380 & 'icntl(24) null pivot detection(0=off) =',I10/
2381 & 'icntl(31) discard factors(0=off, else=on) =',I10/
2382 & 'icntl(32) forward elimination during facto(0=off)=',I10/
2383 & 'icntl(33) compute determinant (0=off) =',I10/
2384 & 'icntl(35) block low rank(blr, 0=off >0=on) =',I10/
2385 & 'icntl(36) blr variant =',I10)
2386 993 FORMAT (
2387 & 'icntl(14) percent of memory increase =',I10)
2388 995 FORMAT (
2389 & 'icntl(9) solve a x=b(1) or a''x = b(else) =',I10/
2390 & 'icntl(10) max steps iterative refinement =',I10/
2391 & 'icntl(11) error analysis(1=all,2=some,else=off) =',I10/
2392 & 'icntl(20) den.(0)/sparse(1,2,3)/dist.(10,11) rhs =',I10/
2393 & 'icntl(21) gathered(0) or distributed(1) solution =',I10)
2394 998 FORMAT (
2395 & ' Size of schur matrix(size_schur) =',I10)
max
#define max(a, b)
Definition macros.h:21

◆ dmumps_print_keep()

subroutine dmumps_print_keep ( type (dmumps_struc), intent(in), target id,
integer lp )

Definition at line 2398 of file dmumps_driver.F.

2399 USE dmumps_struc_def
2400*
2401* ==========
2402* Parameters
2403* ==========
2404 TYPE (DMUMPS_STRUC), TARGET, INTENT(IN) :: id
2405 INTEGER ::LP
2406** Local Variables
2407 INTEGER, POINTER :: JOB
2408 INTEGER,DIMENSION(:),POINTER::ICNTL, KEEP
2409 INTEGER MASTER
2410 parameter( master = 0 )
2411 IF (lp.LE.0) RETURN
2412 job=>id%JOB
2413 icntl=>id%ICNTL
2414 keep=>id%KEEP
2415 IF (id%MYID.EQ.master) THEN
2416 SELECT CASE (job)
2417 CASE(1);
2418 WRITE (lp,980)
2419 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2420 WRITE (lp,991) keep(55),keep(23),icntl(7),keep(95),
2421 & icntl(13),keep(54),keep(60),icntl(22)
2422 IF ((keep(23).EQ.5).OR.(keep(23).EQ.6))THEN
2423 WRITE (lp,992) keep(52)
2424 ENDIF
2425 WRITE (lp,993) keep(12)
2426 CASE(2);
2427 WRITE (lp,980)
2428 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2429 IF (keep(23).EQ.0)THEN
2430 WRITE (lp,992) keep(52)
2431 ENDIF
2432 WRITE (lp,993) keep(12)
2433 CASE(3);
2434 WRITE (lp,980)
2435 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2436 WRITE (lp,995)
2437 & icntl(9),icntl(10),icntl(11),icntl(20),icntl(21)
2438 CASE(4);
2439 WRITE (lp,980)
2440 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2441 IF (keep(23).NE.0)THEN
2442 WRITE (lp,992) keep(52)
2443 ENDIF
2444 WRITE (lp,991) keep(55),keep(23),icntl(7),keep(95),
2445 & icntl(13),keep(54),keep(60),icntl(22)
2446 WRITE (lp,995)
2447 & icntl(9),icntl(10),icntl(11),icntl(20),icntl(21)
2448 WRITE (lp,993) keep(12)
2449 CASE(5);
2450 WRITE (lp,980)
2451 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2452 WRITE (lp,991) keep(55),keep(23),icntl(7),keep(95),
2453 & icntl(13),keep(54),keep(60),icntl(22)
2454 IF ((keep(23).EQ.5).OR.(keep(23).EQ.6)
2455 & .OR. (keep(23).EQ.7)) THEN
2456 WRITE (lp,992) keep(52)
2457 ENDIF
2458 IF (keep(23).EQ.0)THEN
2459 WRITE (lp,992) keep(52)
2460 ENDIF
2461 WRITE (lp,993) keep(12)
2462 CASE(6);
2463 WRITE (lp,980)
2464 WRITE (lp,990) icntl(1),icntl(2),icntl(3),icntl(4)
2465 WRITE (lp,991) keep(55),keep(23),icntl(7),keep(95),
2466 & icntl(13),keep(54),keep(60),icntl(22)
2467 IF ((keep(23).EQ.5).OR.(keep(23).EQ.6)
2468 & .OR. (keep(23).EQ.7)) THEN
2469 WRITE (lp,992) keep(52)
2470 ENDIF
2471 IF (keep(23).EQ.0)THEN
2472 WRITE (lp,992) keep(52)
2473 ENDIF
2474 WRITE (lp,995)
2475 & icntl(9),icntl(10),icntl(11),keep(248),icntl(21)
2476 WRITE (lp,993) keep(12)
2477 END SELECT
2478 ENDIF
2479 980 FORMAT (/'******INTERNAL VALUE OF PARAMETERS (ICNTL/KEEP)****'/)
2480 990 FORMAT (
2481 & 'ICNTL(1) Output stream for error messages =',i10/
2482 & 'ICNTL(2) Output stream for diagnostic messages =',i10/
2483 & 'ICNTL(3) Output stream for global information =',i10/
2484 & 'ICNTL(4) Level of printing =',i10)
2485 991 FORMAT (
2486 & 'ICNTL(5) Matrix format ( keep(55) ) =',i10/
2487 & 'ICNTL(6) Maximum transversal ( keep(23) ) =',i10/
2488 & 'ICNTL(7) Ordering =',i10/
2489 & 'ICNTL(12) LDLT ordering strat ( keep(95) ) =',i10/
2490 & 'ICNTL(13) Parallel root (0=on, 1=off) =',i10/
2491 & 'ICNTL(18) Distributed matrix ( keep(54) ) =',i10/
2492 & 'ICNTL(19) Schur option ( keep(60) 0=off,else=on ) =',i10/
2493 & 'ICNTL(22) Out-of-core option (0=Off, >0=ON) =',i10)
2494 992 FORMAT (
2495 & 'ICNTL(8) Scaling strategy ( keep(52) ) =',i10)
2496 993 FORMAT (
2497 & 'ICNTL(14) Percent of memory increase ( keep(12) ) =',i10)
2498 995 FORMAT (
2499 & 'ICNTL(9) Solve A x=b (1) or A''x = b (else) =',i10/
2500 & 'ICNTL(10) Max steps iterative refinement =',i10/
2501 & 'ICNTL(11) Error analysis ( 0= off, else=on) =',i10/
2502 & 'ICNTL(20) Den.(0)/sparse(1,2,3)/dist.(10,11) RHS =',i10/
2503 & 'ICNTL(21) Gathered (0) or distributed(1) solution =',i10)

◆ dmumps_set_infog()

subroutine dmumps_set_infog ( integer, dimension(80) info,
integer, dimension(size_infog) infog,
integer comm,
integer myid )

Definition at line 2168 of file dmumps_driver.F.

2169 IMPLICIT NONE
2170 include 'mpif.h'
2171C
2172C Purpose:
2173C =======
2174C
2175C If one proc has INFO(1).lt.0 and INFO(1) .ne. -1,
2176C puts INFO(1:2) of this proc on all procs in INFOG
2177C
2178C Arguments:
2179C =========
2180C
2181 INTEGER, PARAMETER :: SIZE_INFOG = 80
2182 INTEGER :: INFO(80)
2183 INTEGER :: INFOG(SIZE_INFOG) ! INFOG(80)
2184 INTEGER :: COMM, MYID
2185C
2186C Local variables
2187C ===============
2188C
2189#if defined(WORKAROUNDINTELILP64MPI2INTEGER)
2190 INTEGER(4) :: TMP1(2),TMP(2)
2191#else
2192 INTEGER :: TMP1(2),TMP(2)
2193#endif
2194 INTEGER ROOT, IERR
2195 INTEGER MASTER
2196 parameter(master=0)
2197C
2198C
2199 IF ( info(1) .ge. 0 ) THEN
2200C
2201C This can only happen if the phase was successful
2202C on all procs. If one proc failed, then all other
2203C procs would have INFO(1)=-1.
2204C
2205 infog(1) = info(1)
2206 infog(2) = info(2)
2207 ELSE
2208C ---------------------
2209C Find who has smallest
2210C error code INFO(1)
2211C ---------------------
2212 infog(1) = info(1)
2213C INFOG(2) = MYID
2214 tmp1(1) = info(1)
2215 tmp1(2) = myid
2216 CALL mpi_allreduce(tmp1,tmp,1,mpi_2integer,
2217 & mpi_minloc,comm,ierr )
2218 infog(2) = info(2)
2219 root = tmp(2)
2220 CALL mpi_bcast( infog(1), 1, mpi_integer, root, comm, ierr )
2221 CALL mpi_bcast( infog(2), 1, mpi_integer, root, comm, ierr )
2222 END IF
2223C
2224C Make INFOG available on all procs:
2225C
2226 CALL mpi_bcast(infog(3), size_infog-2, mpi_integer,
2227 & master, comm, ierr )
2228 RETURN

◆ dmumps_set_k221()

subroutine dmumps_set_k221 ( type (dmumps_struc) id)

Definition at line 2562 of file dmumps_driver.F.

2563 USE dmumps_struc_def
2564 IMPLICIT NONE
2565C
2566C Purpose:
2567C =======
2568C
2569C Sets KEEP(221) on master.
2570C Constraint: must be called before DMUMPS_CHECK_REDRHS.
2571C Can be called at factorization or solve phase
2572C
2573 TYPE (DMUMPS_STRUC) :: id
2574 INTEGER MASTER
2575 parameter( master = 0 )
2576 IF (id%MYID.EQ.master) THEN
2577 id%KEEP(221)=id%ICNTL(26)
2578 IF (id%KEEP(221).ne.0 .and. id%KEEP(221) .NE.1
2579 & .AND.id%KEEP(221).ne.2) id%KEEP(221)=0
2580 ENDIF
2581 RETURN