QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
dqrm_factorize.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 
44 subroutine dqrm_factorize(qrm_mat, transp)
45 
46  use dqrm_spmat_mod
47  use qrm_error_mod
49  use dqrm_fdata_mod
50  use dqrm_utils_mod
51  use qrm_string_mod
52  use qrm_common_mod
53  implicit none
54 
55  type(dqrm_spmat_type), target :: qrm_mat
56  character, optional, intent(in) :: transp
57 
58  integer :: i, totnnz, info, h
59  real(kind(1.d0)) :: t1, t2
60  type(dqrm_front_type), pointer :: front
61  character :: itransp
62  integer, pointer :: tmp(:)
63  ! error management
64  integer :: err_act
65  character(len=*), parameter :: name='qrm_factorize'
66 
67  call qrm_err_act_save(err_act)
68 
69  __qrm_prnt_dbg('("Entering the factorization driver")')
70 
71  ! immediately check if the analysis was done. Otherwise push an error and return
72  if(.not. qrm_mat%adata%ok) then
73  call qrm_err_push(13, 'qrm_factorize')
74  goto 9999
75  end if
76 
77  call dqrm_check_spmat(qrm_mat, qrm_factorize_)
78  __qrm_check_ret(name,'qrm_check_spmat',9999)
79 
80 
81  if(present(transp)) then
82  itransp = transp
83  else
84  itransp = 'n'
85  end if
86 
87  ! in case transp=='t' switch temporarily the row and column indices as well as m and n
88  if(qrm_str_tolower(itransp) .eq. 't') then
89  tmp => qrm_mat%irn
90  qrm_mat%irn => qrm_mat%jcn
91  qrm_mat%jcn => tmp
92  i = qrm_mat%m
93  qrm_mat%m = qrm_mat%n
94  qrm_mat%n = i
95 #if defined(zprec) || defined(cprec)
96  qrm_mat%val = conjg(qrm_mat%val)
97 #endif
98  end if
99 
100  ! initialize the data for the facto
101  call dqrm_factorization_init(qrm_mat)
102  __qrm_check_ret(name,'qrm_factorization_init',9998)
103 
104  !$ call omp_set_num_threads(1)
105 
106  call dqrm_factorization_core(qrm_mat)
107  __qrm_check_ret(name,'qrm_factorization_core',9998)
108 
109 9998 continue
110 
111  qrm_mat%gstats(qrm_nnz_r_) = 0
112  qrm_mat%gstats(qrm_nnz_h_) = 0
113  do i=1, qrm_mat%adata%nnodes
114  qrm_mat%gstats(qrm_nnz_r_) = qrm_mat%gstats(qrm_nnz_r_) + &
115  & qrm_mat%fdata%front_list(i)%rsize
116  end do
117 
118  if(qrm_mat%icntl(qrm_keeph_) .eq. qrm_yes_) then
119  do i=1, qrm_mat%adata%nnodes
120  qrm_mat%gstats(qrm_nnz_h_) = qrm_mat%gstats(qrm_nnz_h_) + &
121  & qrm_mat%fdata%front_list(i)%hsize
122  end do
123  end if
124 
125 
126 
127  ! in case transp=='t' switch temporarily the row and column indices as well as m and n
128  if(qrm_str_tolower(itransp) .eq. 't') then
129  tmp => qrm_mat%irn
130  qrm_mat%irn => qrm_mat%jcn
131  qrm_mat%jcn => tmp
132  i = qrm_mat%m
133  qrm_mat%m = qrm_mat%n
134  qrm_mat%n = i
135 #if defined(zprec) || defined(cprec)
136  qrm_mat%val = conjg(qrm_mat%val)
137 #endif
138  end if
139 
140  call qrm_err_get(info)
141  if(info .gt. 0) goto 9999
142 
143  ! the factorization was succesfully performed
144  qrm_mat%fdata%ok = .true.
145 
146  call qrm_err_act_restore(err_act)
147  return
148 
149 9999 continue ! error management
150  call qrm_err_act_restore(err_act)
151  if(err_act .eq. qrm_abort_) then
152  call qrm_err_check()
153  end if
154 
155  return
156 
157 end subroutine dqrm_factorize
158 
This module contains all the error management routines and data.
This module contains the interfaces of all non-typed routines.
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...
This module contains the definition of the basic sparse matrix type and of the associated methods...
This module contains generic interfaces for a number of auxiliary tools.
subroutine dqrm_factorization_core(qrm_mat)
This is the main factorization routine. It performs the factorization of all the fronts that have bee...
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
This module contains the definition of all the data related to the factorization phase.
This module contains all the generic interfaces for the typed routines in the factorization phase...
This type defines the data structure used to store a matrix.
subroutine dqrm_check_spmat(qrm_spmat, op)
Check the compatibility and correctness of icntl and rcntl parameters.
subroutine dqrm_factorization_init(qrm_mat)
This subroutine initializes the data structures needed for the actual factorization.
This type defines a data structure containing all the data related to a front.
This module contains various string handling routines.
int i
Definition: secs.c:40
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
subroutine dqrm_factorize(qrm_mat, transp)
This routine is the main factorization driver.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.