QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_error_mod.F90
Go to the documentation of this file.
1 !! ##############################################################################################
2 !!
3 !! Copyright 2012 CNRS, INPT
4 !!
5 !! This file is part of qr_mumps.
6 !!
7 !! qr_mumps is free software: you can redistribute it and/or modify
8 !! it under the terms of the GNU Lesser General Public License as
9 !! published by the Free Software Foundation, either version 3 of
10 !! the License, or (at your option) any later version.
11 !!
12 !! qr_mumps is distributed in the hope that it will be useful,
13 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
14 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 !! GNU Lesser General Public License for more details.
16 !!
17 !! You can find a copy of the GNU Lesser General Public License
18 !! in the qr_mumps/doc directory.
19 !!
20 !! ##############################################################################################
21 
22 
23 !! ##############################################################################################
33 
34 
35 #include "qrm_common.h"
36 
39 
90 
92 
95 
97  integer :: code
98 
99  character(len=30) :: sub=''
100 
102  integer, dimension(5) :: i_err_data=0
103 
104  character(len=40) :: a_err_data=''
105 
106  type(qrm_err_type), pointer :: next=>null()
107  end type qrm_err_type
108 
112  type(qrm_err_type), pointer :: top=>null()
113 
114  integer :: nelem=0
115  end type qrm_err_stack_type
116 
118  type(qrm_err_stack_type), save :: qrm_err_stack ! the error messages stack
119 
120  integer, parameter :: qrm_abort_=0, qrm_return_=1
121 
122  integer :: qrm_err_act=qrm_abort_
123 
124 contains
128  subroutine qrm_err_act_save(err_act)
129  integer :: err_act
130  err_act = qrm_err_act
131  qrm_err_act = qrm_return_
132  return
133  end subroutine qrm_err_act_save
134 
135 
138  subroutine qrm_err_act_set(err_act)
139  integer :: err_act
140  if((err_act .ne. qrm_abort_) .and. &
141  & (err_act .ne. qrm_return_)) then
142  call qrm_err_push(26, 'qrm_err_act_save',ied=(/err_act,0,0,0,0/))
143  if(qrm_err_act .eq. qrm_abort_) call qrm_err_check()
144  else
145  qrm_err_act = err_act
146  end if
147 
148  return
149  end subroutine qrm_err_act_set
150 
151 
152 
155  subroutine qrm_err_act_restore(err_act)
156  integer :: err_act
157  qrm_err_act = err_act
158  return
159  end subroutine qrm_err_act_restore
160 
170  subroutine qrm_err_push(code, sub, ied, aed)
171 
172  implicit none
173 
174  integer :: code
175  character(len=*), optional :: sub
176  integer, dimension(5), optional :: ied
177  character(len=*), optional :: aed
178 
179  type(qrm_err_type), pointer :: new_err
180 
181  !$omp critical (error)
182  allocate(new_err)
183 
184  new_err%code = code
185  if(present(sub)) new_err%sub = sub
186  if(present(ied)) new_err%i_err_data = ied
187  if(present(aed)) new_err%a_err_data = aed
188 
189  new_err%next => qrm_err_stack%top
190  qrm_err_stack%top => new_err
191  qrm_err_stack%nelem = qrm_err_stack%nelem+1
192  !$omp end critical (error)
193 
194  return
195 
196  end subroutine qrm_err_push
197 
198 
209  subroutine qrm_err_raise(code, sub, ied, aed)
210  implicit none
211 
212  integer :: code
213  character(len=*), optional :: sub
214  integer, dimension(5), optional :: ied
215  character(len=*), optional :: aed
216 
217  call qrm_err_push(code, sub, ied, aed)
218  call qrm_err_check()
219 
220  return
221 
222  end subroutine qrm_err_raise
223 
224 
227  subroutine qrm_err_get(info)
228 
229  implicit none
230 
231  integer :: info
232  type(qrm_err_type), pointer :: curr
233  integer :: i
234  ! local check
235 
236  info = 0
237  curr => qrm_err_stack%top
238 
239  do i=1, qrm_err_stack%nelem-1
240  curr => curr%next
241  end do
242 
243  if(associated(curr)) info = curr%code
244 
245  return
246  end subroutine qrm_err_get
247 
248 
252  subroutine qrm_err_check( )
253  use qrm_const_mod
254  implicit none
255 
256  ! local check
257  if(qrm_err_stack%nelem .gt. 0) then
258  ! errors were previously detected
259  __qrm_prnt_err('(" ")')
260  __qrm_prnt_err('("==================== Error detected in qr_mumps ====================")')
261  call qrm_flush_err_stack(.true.)
262  __qrm_prnt_err('("====================================================================")')
263  stop
264  end if
265 
266  return
267 
268  end subroutine qrm_err_check
269 
270 
274  subroutine qrm_flush_err_stack(prnt)
275  implicit none
276  logical, optional :: prnt
277 
278  type(qrm_err_type), pointer :: curr, tmp
279  logical :: iprnt
280 
281  if(present(prnt)) then
282  iprnt = prnt
283  else
284  iprnt = .true.
285  end if
286 
287  curr => qrm_err_stack%top
288 
289  do while (associated(curr))
290  if(iprnt) call qrm_process_msg(curr)
291  tmp => curr
292  curr => curr%next
293  deallocate(tmp)
294  qrm_err_stack%nelem = qrm_err_stack%nelem-1
295  end do
296 
297  qrm_err_stack%top => null()
298 
299  return
300 
301  end subroutine qrm_flush_err_stack
302 
303 
309  subroutine qrm_process_msg(msg)
310  use qrm_const_mod
311  implicit none
312 
313  type(qrm_err_type), pointer :: msg
314 
315  __qrm_prnt_err('("Error in subroutine ",a30 " :")')msg%sub
316  select case(msg%code)
317  case(0)
318  __qrm_prnt_err('("Generic error")')
319  case(1)
320  __qrm_prnt_err('("Sparse matrix format ",a3," is not (yet) supported.")')msg%a_err_data(1:3)
321  case(2)
322  __qrm_prnt_err('("Symmetric matrices are not supported.")')
323  case(3)
324  __qrm_prnt_err('("qrm_spmat%cntl is not associated/valid.")')
325  case(4)
326  __qrm_prnt_err('("Trying to allocate an already allocated array.")')
327  case(5)
328  __qrm_prnt_err('("Memory allocation problem. Size required: ",i30)')msg%i_err_data(1)
329  case(6)
330  __qrm_prnt_err('("Trying to deallocate an unallocated array.")')
331  case(7)
332  __qrm_prnt_err('("Memory deallocation problem. ",i30)')msg%i_err_data(1)
333  case(8)
334  __qrm_prnt_err('("Input column permutation not provided/valid")')
335  case(9)
336  __qrm_prnt_err('("Requested ordering method unknown: ",i3)')msg%i_err_data(1)
337  case(10)
338  __qrm_prnt_err('("Insufficient size for array: ",a20)')msg%a_err_data(1:20)
339  case(11)
340  __qrm_prnt_err('("Error in lapack routine: ",i3)')msg%i_err_data(1)
341  case(12)
342  __qrm_prnt_err('("No more memory available")')
343  case(13)
344  __qrm_prnt_err('("The analysis must be done before the factorization")')
345  case(14)
346  __qrm_prnt_err('("The factorization must be done before the solve")')
347  case(15)
348  __qrm_prnt_err('("This type of norm is not implemented.")')
349  case(16)
350  __qrm_prnt_err('("Requested ordering method not available: ",a20)')msg%a_err_data
351  case(17)
352  __qrm_prnt_err('("Error from call to subroutine ",a30)')msg%a_err_data
353  case(18)
354  __qrm_prnt_err('("COLAMD error ")')
355  case(19)
356  __qrm_prnt_err('("SCOTCH error ")')
357  case(20)
358  __qrm_prnt_err('("Factorization error ")')
359  case(21)
360  __qrm_prnt_err('("Apply error ")')
361  case(22)
362  __qrm_prnt_err('("Solve error ")')
363  case(23)
364  __qrm_prnt_err('("Incorrect argument to qrm_set/qrm_get ",a30)')msg%a_err_data
365  case(25)
366  __qrm_prnt_err('("Problem opening file ",a30)')msg%a_err_data
367  case(26)
368  __qrm_prnt_err('("Unknown error action ",i10)')msg%i_err_data(1)
369  case(27)
370  __qrm_prnt_err('("Incompatible values in qrm_spmat%icntl ",i2,2x,i2)')msg%i_err_data(1:2)
371  case(28)
372  __qrm_prnt_err('("Incorrect value for qrm_nb_/qrm_ib_ ",i2,2x,i2)')msg%i_err_data(1)
373  case(29)
374  __qrm_prnt_err('("Incorrect value for qrm_spmat%m/n/nz ",i6,2x,i6,2x,i16)')msg%i_err_data(1:3)
375  case(30)
376  __qrm_prnt_err('("qrm_apply cannot be called if the H matrix is discarded.")')
377  case default
378  __qrm_prnt_err('("Unknown error code",i4)')msg%code
379  end select
380 
381  __qrm_prnt_err('(" ")')
382  return
383 
384  end subroutine qrm_process_msg
385 
386 
387 
388 end module qrm_error_mod
This module contains all the error management routines and data.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine qrm_err_get(info)
This subroutine return the code of the first error on the stack or zero if the stack is empty...
subroutine qrm_err_act_set(err_act)
Sets the default error action.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
This type is to represent the errors stack.
subroutine qrm_process_msg(msg)
This routine prints out a message on the error unit.
subroutine qrm_flush_err_stack(prnt)
This subroutine flushes the errors stack optionally printing all the messages on the eunit output uni...
int i
Definition: secs.c:40
subroutine qrm_err_raise(code, sub, ied, aed)
Pushes an error on the stack and the flushes the stack itself. Basically does err_push and err_check ...
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
This is the basic type for error message.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.