QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
dqrm_factorization_init.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 
38 
50 subroutine dqrm_factorization_init(qrm_mat)
51 
52  use qrm_mem_mod
53  use dqrm_spmat_mod
54  use dqrm_fdata_mod
55  use dqrm_utils_mod
56  use qrm_common_mod
57  use qrm_error_mod
58  implicit none
59 
60  type(dqrm_spmat_type), target :: qrm_mat
61 
62  integer :: f, nrowsf, nvalsf, i, j, roff, r, c, lrow, itmp, k, m, n
63  integer, allocatable :: rmap(:), rcnt(:), vcnt(:), row_to_frow(:)
64  type(dqrm_front_type), pointer :: front
65  real(kind(1.d0)) :: vtmp
66  ! error management
67  integer :: err_act
68  character(len=*), parameter :: name='qrm_analyse'
69 
70  call qrm_err_act_save(err_act)
71 
72  call dqrm_fdata_destroy(qrm_mat%fdata)
73  __qrm_check_ret(name,'qrm_fdata_destroy',9999)
74  allocate(qrm_mat%fdata%front_list(qrm_mat%adata%nnodes))
75 
76  ! rmap contains a mapping of the rows onto fronts
77  call qrm_aalloc(rmap, qrm_mat%m)
78  call qrm_aalloc(row_to_frow, qrm_mat%m)
79  call qrm_aalloc(rcnt, qrm_mat%m)
80  call qrm_aalloc(vcnt, qrm_mat%adata%nnodes)
81  __qrm_check_ret(name,'qrm_aalloc',9999)
82  rcnt = 0
83  vcnt = 0
84  m = qrm_mat%m
85  n = qrm_mat%n
86 
87  qrm_mat%fdata%nfronts = qrm_mat%adata%nnodes
88 #if defined(_OPENMP)
89  call omp_init_lock(qrm_mat%fdata%lock)
90 #endif
91 
92  ! build the row mapping
93  ! TODO: stair(i) may be set to be the first row of front i instead
94  ! of the last
95  roff = 1
96  do f = 1, qrm_mat%adata%nnodes
97  ! for all the rows belonging to this front
98  do i = roff, qrm_mat%adata%stair(f)
99  rmap(qrm_mat%adata%rperm(i)) = f
100  row_to_frow(qrm_mat%adata%rperm(i)) = i-roff+1
101  end do
102  roff = qrm_mat%adata%stair(f)+1
103  end do
104 
105  ! first pass to count
106  do k=1, qrm_mat%nz
107  j = qrm_mat%jcn(k)
108  i = qrm_mat%irn(k)
109  if((j.le.0) .or. (j.gt.n) .or. (i.le.0) .or. (i.gt.m) ) cycle
110  f = rmap(i)
111  ! count the number of values per row
112  rcnt(i) = rcnt(i)+1
113  ! count the number of nonzeros from the original matrix in front f
114  vcnt(f) = vcnt(f)+1
115  end do
116 
117  ! pass to allocate
118  roff = 0
119  qrm_mat%fdata%front_list(:)%status=0
120  do f = 1, qrm_mat%adata%nnodes
121  front => qrm_mat%fdata%front_list(f)
122  nrowsf = qrm_mat%adata%stair(f) - roff
123  front%anrows = nrowsf
124  front%num = f
125  do i=qrm_mat%adata%childptr(f), qrm_mat%adata%childptr(f+1)-1
126  c = qrm_mat%adata%child(i)
127  if(qrm_mat%adata%small(c) .eq. 0) front%status = front%status-1
128  end do
129 #if defined(_OPENMP)
130  call omp_init_lock(front%lock)
131 #endif
132  call qrm_aalloc(front%aiptr, max(nrowsf+1,2))
133  call qrm_aalloc(front%ajcn, vcnt(f))
134  call qrm_aalloc(front%aval, vcnt(f))
135  __qrm_check_ret(name,'qrm_aalloc',9999)
136 
137  front%aiptr(1:2) = 1
138  do i = 1, nrowsf-1
139  front%aiptr(i+2) = front%aiptr(i+1)+rcnt(qrm_mat%adata%rperm(roff+i))
140  end do
141 
142  roff = qrm_mat%adata%stair(f)
143  end do
144 
145 
146 
147  ! pass to fill
148  do k=1, qrm_mat%nz
149  j = qrm_mat%jcn(k)
150  i = qrm_mat%irn(k)
151  if((j.le.0) .or. (j.gt.n) .or. (i.le.0) .or. (i.gt.m) ) cycle
152  r = i
153  c = j
154  f = rmap(r)
155  lrow = row_to_frow(r)
156  front => qrm_mat%fdata%front_list(f)
157  front%ajcn(front%aiptr(lrow+1)) = c
158  front%aval(front%aiptr(lrow+1)) = qrm_mat%val(k)
159  front%aiptr(lrow+1) = front%aiptr(lrow+1)+1
160  end do
161 
162  qrm_mat%fdata%done = 0
163 
164  if(qrm_mat%adata%ncsing .gt. 0) then
165  qrm_mat%fdata%done = 1
166  ! the first front may contain all the singletons. For this
167  ! reason, each diagonal element must be the first of its
168  ! row. This is done to make the solve easier and speed it -up.
169  ! we just sweep over all the elements of the front and swap the
170  ! diagonal element with the first along its row. It may be
171  ! optimized TODO
172  front => qrm_mat%fdata%front_list(1)
173  do i=1, front%anrows
174  do j=front%aiptr(i), front%aiptr(i+1)-1
175  c = front%ajcn(j)
176  if(qrm_mat%adata%cperm(i) .eq. c) then
177  ! this is a diagonal element; swap
178  vtmp = front%aval(front%aiptr(i))
179  front%aval(front%aiptr(i)) = front%aval(j)
180  front%aval(j) = vtmp
181  itmp = front%ajcn(front%aiptr(i))
182  front%ajcn(front%aiptr(i)) = c
183  front%ajcn(j) = itmp
184  end if
185  end do
186  end do
187  end if
188 
189 
190 
191  call qrm_adealloc(rmap)
192  call qrm_adealloc(row_to_frow)
193  call qrm_adealloc(rcnt)
194  call qrm_adealloc(vcnt)
195  __qrm_check_ret(name,'qrm_adealloc',9999)
196 
197  call qrm_err_act_restore(err_act)
198  return
199 
200 9999 continue ! error management
201  call qrm_err_act_restore(err_act)
202  if(err_act .eq. qrm_abort_) then
203  call qrm_err_check()
204  end if
205 
206  return
207 
208 end subroutine dqrm_factorization_init
This module contains all the error management routines and data.
Generic interface for the qrm_adealloc_i, qrm_adealloc_2i, qrm_adealloc_s, qrm_adealloc_2s, qrm_adealloc_3s, qrm_adealloc_d, qrm_adealloc_2d, qrm_adealloc_3d, qrm_adealloc_c, qrm_adealloc_2c, qrm_adealloc_3c, qrm_adealloc_z, qrm_adealloc_2z, qrm_adealloc_3z, routines.
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.
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 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.
Generic interface for the qrm_aalloc_i, qrm_aalloc_2i, qrm_aalloc_s, qrm_aalloc_2s, qrm_aalloc_3s, qrm_aalloc_d, qrm_aalloc_2d, qrm_aalloc_3d, qrm_aalloc_c, qrm_aalloc_2c, qrm_aalloc_3c, qrm_aalloc_z, qrm_aalloc_2z, qrm_aalloc_3z, routines.
Definition: qrm_mem_mod.F90:78
This type defines the data structure used to store a matrix.
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 implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
int i
Definition: secs.c:40
subroutine dqrm_fdata_destroy(qrm_fdata)
Destroys a dqrm_fdata_type instance.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.