35 #include "qrm_common.h"
58 real(kind(1.d0)),
intent(inout) :: b(:,:)
59 real(kind(1.d0)),
intent(out) :: x(:,:)
61 integer :: nth, thn, info, f, dones
66 type(dqrm_front_type
),
pointer :: front
67 integer,
allocatable :: status(:)
68 type(qrm_adata_type
),
pointer :: adata
69 type(dqrm_fdata_type
),
pointer :: fdata
72 integer(kind=omp_lock_kind),
allocatable :: locks(:)
73 integer(kind=omp_lock_kind) :: dlock
78 character(len=*),
parameter :: name=
'qrm_solve_r'
82 __qrm_prnt_dbg(
'("Solving for R")')
85 adata => qrm_mat%adata
86 fdata => qrm_mat%fdata
95 __qrm_check_ret(name,
'qrm_aalloc',9999)
98 do f = 1, adata%nnodes
99 status(f) = qrm_ready_
102 if(p .eq. 0 .and. (adata%rc(f).ge.0)) call
qrm_queue_push(ready_q, f)
107 if(adata%ncsing .gt. 0)
then
121 #if defined (_OPENMP)
122 nth = omp_get_num_threads()
123 thn = omp_get_thread_num()
135 if(qrm_err_stack%nelem .gt. 0) goto 9998
144 if(.not. got_task) cycle taskloop
167 if(qrm_err_stack%nelem .gt. 0)
then
173 if (qrm_mat%adata%ncsing .gt. 0)
then
175 __qrm_check_ret(name,
'qrm_solve_sing_front',9999)
183 if(err_act .eq. qrm_abort_)
then
197 type(dqrm_front_type
),
pointer :: front
204 thn=omp_get_thread_num()
217 front => fdata%front_list(f)
220 if(.not. omp_test_lock(locks(f))) cycle
222 if(status(f) .eq. qrm_ready_)
then
228 status(f) = qrm_busy_
232 call omp_unset_lock(locks(f))
260 if(dones .eq. fdata%nfronts)
then
281 type(dqrm_front_type
),
pointer :: front
282 integer :: f, p, c, info
287 front => qrm_mat%fdata%front_list(task%front)
291 status(task%front) = qrm_done_
295 call omp_set_lock( dlock )
299 call omp_unset_lock( dlock )
304 do p = adata%childptr(front%num), adata%childptr(front%num+1)-1
306 if(adata%small(c) .eq. 1)
then
308 if(info .ne. 0) goto 9997
329 integer :: fnum, info
331 type(dqrm_front_type
),
pointer :: front
332 integer :: node, c, acc, thn, p
339 thn = omp_get_thread_num()
352 front => qrm_mat%fdata%front_list(node)
355 status(node) = qrm_done_
363 do p = adata%childptr(front%num), adata%childptr(front%num+1)-1
393 type(dqrm_front_type
) :: front
396 integer :: k, m, pk, jp, n,
i, thn, cnt, j, pv1, pv2, r
399 real(kind(1.d0)),
allocatable :: in_x(:,:)
402 character(len=*),
parameter :: name=
'front_r'
405 if (min(front%m, front%n) .le.0) goto 9999
406 if (front%npiv .le. 0) goto 9999
409 thn = omp_get_thread_num()
416 call
qrm_aalloc(in_x, front%n,
size(x,2), info)
417 __qrm_check_ret(name,
'qrm_aalloc',9999)
419 in_x(1:front%npiv,:) = b(front%rows(1:front%npiv),:)
420 in_x(front%npiv+1:front%n,:) = x(front%cols(front%npiv+1:front%n),:)
453 outer:
do jp = front%npiv - mod(front%npiv, front%nb)+1, 1, -front%nb
454 pk = min(front%nb, front%npiv-jp+1)
458 inner:
do j = jp+pk-mod(pk,front%ib), jp, -front%ib
459 m = min(front%ib, jp+pk - j)
462 cnt = cnt-m*(front%n-j+1)
464 if(k.gt.0) call dgemm(
'n',
'n', m, n, k, -1.d0, front%r(cnt+m*m), m, &
465 & in_x(j+m,1), front%n, 1.d0, in_x(j,1), front%n)
466 call dtrsm(
'l',
'u',
'n',
'n', m, n, 1.d0, front%r(cnt), &
467 & m, in_x(j,1), front%n)
472 x(front%cols(1:front%npiv),:) = in_x(1:front%npiv,:)
475 __qrm_check_ret(name,
'qrm_adealloc',9999)
This module contains all the error management routines and data.
subroutine qrm_init_task_queue(h)
Inititalizes a set of queues attached to a family of threads referenced through the handle h...
integer function qrm_queue_next(q, n)
Returns the element that follows n in the queue q. Very useful for sweeping through a queue...
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_clean_task_queue(h)
Destroyes a set of queues.
subroutine check_solver_over()
This type defines the handle for the queues attached to a family of threads.
subroutine dqrm_solve_sing_front(qrm_mat, b, x, trans)
This function handles the front containing the singletons during the solve for R or R'...
This module contains all the interfaces for the typed routines in the solve phase.
subroutine solve_r(task, thn)
A data type meant to to define a queue.
logical function qrm_get_task(h, tsk)
Pops a task from a queue. Tasks are always popped from the head of the queue. The return value is ...
subroutine qrm_queue_rm(q, n)
Removes (without returning it) an element from a queue.
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_queue_push(q, elem)
Pushes an element on a queue.
This module contains all the facilities for front queues.
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...
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.
subroutine dqrm_solve_r(qrm_mat, b, x)
This function solves for R against multiple vectors.
subroutine front_r(front, info)
subroutine fill_queue_r()
This type defines a computational task.
subroutine qrm_par_mem_init()
This routine has to be called at the beginning of a parallel section. Afterwards, each thread will up...
integer function qrm_queue_pop(q)
Pops an element from a queue.
This type defines the data structure used to store a matrix.
integer function qrm_task_queue_card(h)
Returns the number of tasks present on a set of queues referenced by a handle.
subroutine qrm_par_mem_finalize()
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
subroutine do_subtree_r(fnum, info)
This module contains an implementation of some operations on triangular/trapezoidal matrices stored i...
logical function qrm_sched_task(h, tsk, pol, q)
Pushes a task on a queue.
This module contains the definition of a task type that is used for scheduling tasks during the facto...
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
subroutine qrm_queue_free(q)
Frees a queue.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.