QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_adata_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 
37 
39 
40 
41 
43 
45 
53 
60  integer, allocatable :: cperm(:)
61 
63  integer, allocatable :: icperm(:)
64 
66  integer, allocatable :: rperm(:)
67 
70  integer, allocatable :: cp_ptr(:)
71 
73  integer, allocatable :: rc(:)
74 
76  integer, allocatable :: parent(:)
77 
78  integer, allocatable :: child(:)
79 
81  integer, allocatable :: childptr(:)
82 
84  integer, allocatable :: nfrows(:)
85 
97  integer, allocatable :: stair(:)
98 
100  integer, allocatable :: small(:)
101 
102  integer, allocatable :: fcol(:)
103 
105  integer, allocatable :: fcol_ptr(:)
106 
108  integer, allocatable :: leaves(:)
109 
110  integer :: nleaves=0
111 
113  integer :: nnodes=0
114 
115  integer :: ncsing=0
116 
117  integer :: nrsing=0
118 
119  logical :: ok=.false.
120  end type qrm_adata_type
121 
122  contains
123 
124 
128  subroutine qrm_adata_copy(adata_in, adata_out)
129 
130  use qrm_mem_mod
131  use qrm_error_mod
132  implicit none
133 
134  type(qrm_adata_type), intent(in) :: adata_in
135  type(qrm_adata_type), intent(out) :: adata_out
136 
137  ! error management
138  integer :: err_act
139  character(len=*), parameter :: name='qrm_adata_copy'
140 
141  call qrm_err_act_save(err_act)
142 
143  call qrm_aalloc(adata_out%cperm, qrm_asize(adata_in%cperm))
144  call qrm_aalloc(adata_out%icperm, qrm_asize(adata_in%icperm))
145  call qrm_aalloc(adata_out%rperm, qrm_asize(adata_in%rperm))
146  call qrm_aalloc(adata_out%cp_ptr, qrm_asize(adata_in%cp_ptr))
147  call qrm_aalloc(adata_out%rc, qrm_asize(adata_in%rc))
148  call qrm_aalloc(adata_out%parent, qrm_asize(adata_in%parent))
149  call qrm_aalloc(adata_out%child, qrm_asize(adata_in%child))
150  call qrm_aalloc(adata_out%childptr, qrm_asize(adata_in%childptr))
151  call qrm_aalloc(adata_out%nfrows, qrm_asize(adata_in%nfrows))
152  call qrm_aalloc(adata_out%stair, qrm_asize(adata_in%stair))
153  call qrm_aalloc(adata_out%leaves, qrm_asize(adata_in%leaves))
154  call qrm_aalloc(adata_out%fcol, qrm_asize(adata_in%fcol))
155  call qrm_aalloc(adata_out%fcol_ptr, qrm_asize(adata_in%fcol_ptr))
156  call qrm_aalloc(adata_out%small, qrm_asize(adata_in%small))
157  __qrm_check_ret(name,'qrm_alloc',9999)
158 
159  if(allocated(adata_in%cperm))adata_out%cperm = adata_in%cperm
160  if(allocated(adata_in%icperm))adata_out%icperm = adata_in%icperm
161  if(allocated(adata_in%rperm))adata_out%rperm = adata_in%rperm
162  if(allocated(adata_in%cp_ptr))adata_out%cp_ptr = adata_in%cp_ptr
163  if(allocated(adata_in%rc))adata_out%rc = adata_in%rc
164  if(allocated(adata_in%parent))adata_out%parent = adata_in%parent
165  if(allocated(adata_in%child))adata_out%child = adata_in%child
166  if(allocated(adata_in%childptr))adata_out%childptr = adata_in%childptr
167  if(allocated(adata_in%nfrows))adata_out%nfrows = adata_in%nfrows
168  if(allocated(adata_in%stair))adata_out%stair = adata_in%stair
169  if(allocated(adata_in%leaves))adata_out%leaves = adata_in%leaves
170  if(allocated(adata_in%fcol))adata_out%fcol = adata_in%fcol
171  if(allocated(adata_in%fcol_ptr))adata_out%fcol_ptr = adata_in%fcol_ptr
172  if(allocated(adata_in%small))adata_out%small = adata_in%small
173  adata_out%nnodes = adata_in%nnodes
174  adata_out%ncsing = adata_in%ncsing
175  adata_out%nrsing = adata_in%nrsing
176  adata_out%nleaves = adata_in%nleaves
177  adata_out%ok = adata_in%ok
178 
179  call qrm_err_act_restore(err_act)
180  return
181 
182 9999 continue ! error management
183  call qrm_err_act_restore(err_act)
184  if(err_act .eq. qrm_abort_) then
185  call qrm_err_check()
186  end if
187  return
188 
189  end subroutine qrm_adata_copy
190 
191 
195  subroutine qrm_adata_move(adata_in, adata_out)
196 
197  use qrm_mem_mod
198  use qrm_error_mod
199  implicit none
200 
201  type(qrm_adata_type), intent(inout) :: adata_in
202  type(qrm_adata_type), intent(inout) :: adata_out
203 
204  ! error management
205  integer :: err_act
206  character(len=*), parameter :: name='qrm_adata_move'
207 
208  call qrm_err_act_save(err_act)
209 
210  call move_alloc(adata_in%cperm , adata_out%cperm )
211  call move_alloc(adata_in%icperm , adata_out%icperm )
212  call move_alloc(adata_in%rperm , adata_out%rperm )
213  call move_alloc(adata_in%cp_ptr , adata_out%cp_ptr )
214  call move_alloc(adata_in%rc , adata_out%rc )
215  call move_alloc(adata_in%parent , adata_out%parent )
216  call move_alloc(adata_in%child , adata_out%child )
217  call move_alloc(adata_in%childptr, adata_out%childptr )
218  call move_alloc(adata_in%nfrows , adata_out%nfrows )
219  call move_alloc(adata_in%stair , adata_out%stair )
220  call move_alloc(adata_in%leaves , adata_out%leaves )
221  call move_alloc(adata_in%fcol , adata_out%fcol )
222  call move_alloc(adata_in%fcol_ptr, adata_out%fcol_ptr )
223  call move_alloc(adata_in%small , adata_out%small )
224 
225  adata_out%nnodes = adata_in%nnodes
226  adata_out%ncsing = adata_in%ncsing
227  adata_out%nrsing = adata_in%nrsing
228  adata_out%nleaves = adata_in%nleaves
229  adata_out%ok = adata_in%ok
230 
231  call qrm_err_act_restore(err_act)
232  return
233 
234 9999 continue ! error management
235  call qrm_err_act_restore(err_act)
236  if(err_act .eq. qrm_abort_) then
237  call qrm_err_check()
238  end if
239  return
240 
241  end subroutine qrm_adata_move
242 
243 
246  subroutine qrm_adata_destroy(adata)
247 
248  use qrm_mem_mod
249  use qrm_error_mod
250  implicit none
251 
252  type(qrm_adata_type) :: adata
253 
254  ! error management
255  integer :: err_act
256  character(len=*), parameter :: name='qrm_adata_destroy'
257 
258  call qrm_err_act_save(err_act)
259 
260  call qrm_adealloc(adata%cperm)
261  call qrm_adealloc(adata%icperm)
262  call qrm_adealloc(adata%rperm)
263  call qrm_adealloc(adata%cp_ptr)
264  call qrm_adealloc(adata%rc)
265  call qrm_adealloc(adata%parent)
266  call qrm_adealloc(adata%nfrows)
267  call qrm_adealloc(adata%stair)
268  call qrm_adealloc(adata%small)
269  call qrm_adealloc(adata%childptr)
270  call qrm_adealloc(adata%child)
271  call qrm_adealloc(adata%leaves)
272  call qrm_adealloc(adata%fcol)
273  call qrm_adealloc(adata%fcol_ptr)
274  __qrm_check_ret(name,'qrm_dealloc',9999)
275 
276  adata%nleaves = 0
277  adata%nnodes = 0
278  adata%ncsing = 0
279  adata%nrsing = 0
280 
281  adata%ok = .false.
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_adata_destroy
293 
294 
295 end module qrm_adata_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.
subroutine qrm_adata_copy(adata_in, adata_out)
make s a copy of an qrm_adata_type instance
subroutine qrm_adata_move(adata_in, adata_out)
make s a move of an qrm_adata_type instance
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
This module contains the definition of the analysis data type.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine qrm_adata_destroy(adata)
Frees an qrm_adata_type instance.
The main data type for the analysis 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
Generic interface for the qrm_asize_i, qrm_asize_s, qrm_asize_2s, qrm_asize_3s, qrm_asize_d, qrm_asize_2d, qrm_asize_3d, qrm_asize_c, qrm_asize_2c, qrm_asize_3c, qrm_asize_z, qrm_asize_2z, qrm_asize_3z routines.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.