35 #include "qrm_common.h"
54 integer ::
i, j, row1, col1, row2, col2
55 integer,
allocatable :: mark(:)
58 character(len=*),
parameter :: name=
'qrm_ata_graph'
66 ata_graph%iptr(1:2) = 1
69 __qrm_check_ret(name,
'convert/alloc',9999)
75 do i=g_csc%jptr(col1), g_csc%jptr(col1+1)-1
80 do j=g_csr%iptr(row1), g_csr%iptr(row1+1)-1
85 if(col1 .eq. col2) cycle
86 if(mark(col2) .lt. col1)
then
88 ata_graph%iptr(col1+2) = ata_graph%iptr(col1+2)+1
95 ata_graph%iptr(
i) = ata_graph%iptr(
i)+ata_graph%iptr(
i-1)
98 ata_graph%nz = ata_graph%iptr(g_csc%n+2)
100 __qrm_check_ret(name,
'qrm_palloc',9999)
106 do i=g_csc%jptr(col1), g_csc%jptr(col1+1)-1
111 do j=g_csr%iptr(row1), g_csr%iptr(row1+1)-1
116 if(col1 .eq. col2) cycle
117 if(mark(col2) .lt. col1)
then
119 ata_graph%jcn(ata_graph%iptr(col1+1)) = col2
120 ata_graph%iptr(col1+1) = ata_graph%iptr(col1+1)+1
126 ata_graph%n = g_csc%n
127 ata_graph%m = g_csc%n
131 __qrm_check_ret(name,
'destroy/dealloc',9999)
138 if(err_act .eq. qrm_abort_)
then
This module contains all the error management routines and data.
subroutine _qrm_spmat_convert(in_mat, out_mat, fmt, values)
This subroutine converts an input matrix into a different storage format. Optionally the values may b...
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_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
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.
This type defines the data structure used to store a matrix.
This module contains the definition of the basic sparse matrix type and of the associated methods...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Generic interface for the qrm_palloc_i, qrm_palloc_2i, qrm_palloc_s, qrm_palloc_2s, qrm_palloc_d, qrm_palloc_2d, qrm_palloc_c, qrm_palloc_2c, qrm_palloc_z, qrm_palloc_2z, routines.
subroutine _qrm_ata_graph(g_csc, ata_graph)
This subroutine computes the fill reducing ordering using METIS.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.