QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_fdata_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 
41 
42 #if defined(_OPENMP)
43  use omp_lib
44 #endif
45 
46  integer, parameter :: qrm_busy_ = -huge(0)
47  integer, parameter :: qrm_ready_ = 0 ! this can never be changed
48  integer, parameter :: qrm_activable_ = 1
49  integer, parameter :: qrm_active_ = 2
50  integer, parameter :: qrm_factorized_ = 3 ! the front has been completely factorized
51  integer, parameter :: qrm_free_ = 4 ! the front has been freed
52  integer, parameter :: qrm_done_ = huge(0)
53 
57  integer :: num=0
58 
59  integer :: m=0
60 
61  integer :: n=0
62 
63  integer :: npiv=0
64 
65  integer, allocatable :: rows(:)
66 
67  integer, allocatable :: cols(:)
68 
70  integer, allocatable :: aiptr(:)
71 
73  integer, allocatable :: ajcn(:)
74 
77  _qrm_data, allocatable :: aval(:)
79  ! integer contains the number of rows in this front from the original matrix
80  integer :: anrows
81 
84  integer, allocatable :: colmap(:)
85 
90  integer, allocatable :: rowmap(:)
91 
94  integer, allocatable :: stair(:)
95 
97  _qrm_data, allocatable :: front(:,:)
100  _qrm_data, allocatable :: r(:)
103  _qrm_data, allocatable :: h(:)
107  _qrm_data, allocatable :: cb(:)
110  _qrm_data, allocatable :: tau(:)
113  _qrm_data, allocatable :: t(:,:)
116  _qrm_data, allocatable :: b(:,:)
127  integer, allocatable :: ptable(:)
128 
140  integer :: status
141 
143  integer :: owner
144 
145  integer :: lastpnl=0
146 
147  integer :: nb, ib
148 
149  integer :: nc
150 
151  integer :: np
152 
154  integer :: ne
155 
156  integer(kind=8) :: rsize, hsize
157 #if defined(_OPENMP)
158 
159  integer(kind=omp_lock_kind) :: lock
160 #endif
161  end type _qrm_front_type
162 
163 
164 
169  integer :: nfronts=0
170 
172  integer :: done=0
173 
175  type(_qrm_front_type), allocatable :: front_list(:)
176 #if defined(_OPENMP)
177 
178  integer(kind=omp_lock_kind) :: lock
179 #endif
180 
181  logical :: ok=.false.
182  end type _qrm_fdata_type
183 
185  module procedure _qrm_fdata_destroy
186  end interface qrm_fdata_destroy
187 
188 
189 
190 contains
191 
192  ! TODO: add copy!
193 
198  subroutine _qrm_front_destroy(qrm_front)
199  use qrm_mem_mod
200  use qrm_error_mod
201  implicit none
202 
203  type(_qrm_front_type) :: qrm_front
204 
205  ! error management
206  integer :: err_act
207  character(len=*), parameter :: name='qrm_front_destroy'
208 
209  call qrm_err_act_save(err_act)
210 
211  call qrm_adealloc(qrm_front%aiptr)
212  call qrm_adealloc(qrm_front%ajcn)
213  call qrm_adealloc(qrm_front%aval)
214  call qrm_adealloc(qrm_front%front)
215  call qrm_adealloc(qrm_front%r)
216  call qrm_adealloc(qrm_front%h)
217  call qrm_adealloc(qrm_front%rows)
218  call qrm_adealloc(qrm_front%cols)
219  call qrm_adealloc(qrm_front%rowmap)
220  call qrm_adealloc(qrm_front%colmap)
221  call qrm_adealloc(qrm_front%stair)
222  call qrm_adealloc(qrm_front%ptable)
223  call qrm_adealloc(qrm_front%cb)
224  call qrm_adealloc(qrm_front%tau)
225  call qrm_adealloc(qrm_front%t)
226  call qrm_adealloc(qrm_front%b)
227 #if defined(_OPENMP)
228  call omp_destroy_lock(qrm_front%lock)
229 #endif
230  __qrm_check_ret(name,'qrm_adelloc',9999)
231 
232  qrm_front%m = 0
233  qrm_front%n = 0
234 
235  call qrm_err_act_restore(err_act)
236  return
237 
238 9999 continue ! error management
239  call qrm_err_act_restore(err_act)
240  if(err_act .eq. qrm_abort_) then
241  call qrm_err_check()
242  end if
243  return
244 
245  end subroutine _qrm_front_destroy
246 
247 
248 
253  subroutine _qrm_fdata_destroy(qrm_fdata)
254  use qrm_mem_mod
255  use qrm_error_mod
256  implicit none
257 
258  type(_qrm_fdata_type) :: qrm_fdata
259  integer :: i
260  ! error management
261  integer :: err_act
262  character(len=*), parameter :: name='qrm_fdata_destroy'
263 
264  call qrm_err_act_save(err_act)
265 
266 
267  if(allocated(qrm_fdata%front_list)) then
268  do i=1, qrm_fdata%nfronts
269  call _qrm_front_destroy(qrm_fdata%front_list(i))
270  end do
271  deallocate(qrm_fdata%front_list)
272  end if
273  __qrm_check_ret(name,'qrm_front_destroy',9999)
274 
275  qrm_fdata%done = 0
276  qrm_fdata%nfronts = 0
277  qrm_fdata%ok = .false.
278 #if defined(_OPENMP)
279  ! FIXME: how to check whether a lock is initialized or not?
280  ! call omp_destroy_lock(qrm_fdata%lock)
281 #endif
282  call qrm_err_act_restore(err_act)
283  return
284 
285 9999 continue ! error management
286  call qrm_err_act_restore(err_act)
287  if(err_act .eq. qrm_abort_) then
288  call qrm_err_check()
289  end if
290  return
291 
292  end subroutine _qrm_fdata_destroy
293 
294 end module _qrm_fdata_mod
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 definition of all the data related to the factorization phase.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine _qrm_fdata_destroy(qrm_fdata)
Destroys a _qrm_fdata_type instance.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine _qrm_front_destroy(qrm_front)
Frees a qrm_front_type instance.
The data structure meant to store all the results of the factorization phase.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
int i
Definition: secs.c:40
This type defines a data structure containing all the data related to a front.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.