Compare commits

..
Author SHA1 Message Date
Stack-1 3a5e487d66 [ADD] Warm-up repetitions in the AMG comm-scheme driver
Fifth positional argument nwarm (default 1). The first run of a scheme
pays the one-off initialisation of the one-sided transport: ~1600 s
aggregated at 448 ranks, charged entirely to whichever scheme creates
the first window. Warm-up rows are written to the CSV as
warmup_uniform / warmup_sensitivity, so the setup cost stays in the
data while staying out of the averages.
2026-08-13 18:02:57 +02:00
Stack-1 9168d450ce [ADD] Added scripts to collect data 2026-08-11 10:44:53 +02:00
Stack-1 3839ab8baa [ADD] Added first test for communication schemes profiling on additive and multiplicative AMG 2026-08-09 17:51:24 +02:00
Stack-1 7030eaea64 [UPDATE] Removed empty missing file 2026-06-15 13:18:50 +02:00
Stack-1 4031ffb7ba [UPDATE] Drop work= from the vector (_vect) apply chain for communication_v2
Adapt amg4psblas to the communication_v2 PSBLAS interfaces, which removed
the work argument from all vector routines (psb_*_vect, the base prec
apply_vect, map_U2V_v/map_V2U_v).

Remove work from the whole vector apply chain:
- prec apply2_vect/apply1_vect and the amg_*precaply2/1_vect implementations
  (dropped the local work_ buffer);
- mlprec_aply_vect and its inner recursive routines
  (inner_ml_aply, inner_add/mult/k_cycle, inneritkcycle);
- the smoother/solver *_apply_vect implementations and their interface
  declarations, dropping the now-dead aux/ww scratch buffers;
- AS smoother restr_v/prol_v (psb_halo/psb_ovrl on vectors);
- onelev map_rstr_v/map_prol_v (psb_map_U2V/V2U on vectors);
- poly_smoother_bld power-iteration apply_v calls.

Array routines keep work. mumps/slu/sludist/umf vector applies now pass a
zero-size local buffer to the underlying array apply (work optional there).
Regenerated for s/d/c/z; the library builds and the pdegen samples converge.
2026-06-15 12:45:11 +02:00
sfilippone a3e1be46ee Fix matrix generation 2026-05-05 16:45:06 +02:00
sfilippone 3343b039e6 Fix sample matrix generators 2026-05-05 15:52:03 +02:00
sfilippone 5c055170e7 Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-05-05 15:34:47 +02:00
sfilippone 0814492adc Mods jac solver 2026-05-05 15:33:45 +02:00
sfilippone 1ae3cc135f Improve error handling for prec%free 2026-05-04 20:50:34 +02:00
sfilippone 246992cb65 Add base matrix info to prec%descr 2026-04-27 13:04:50 +02:00
sfilippone 01cc7ada88 Adjust strategy for stopping on aggregation ratio 2026-04-27 13:04:04 +02:00
sfilippone 693eab66cb Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-04-16 14:59:33 +02:00
sfilippone 1a2ec161d7 Improve configry for MUMPS 2026-04-16 14:59:11 +02:00
sfilippone 4642c857d1 Improve MUMPS solver build 2026-04-16 14:58:59 +02:00
sfilippone c8d065fa55 Fix double allocation in DDIAG%BLD 2026-04-15 10:19:50 +02:00
sfilippone 0e1d7de857 Update VERSION 2026-04-10 13:37:51 +02:00
Salvatore Filippone a22725b787 CMakeLists interim version 2026-04-09 14:37:15 +02:00
sfilippone 74a9ca90cb Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-04-07 16:24:04 +02:00
sfilippone 0d38955a2d Fix generation of amg_config.h with CMAKE 2026-04-07 16:23:32 +02:00
fdurastante 9be4a3f1be Updated MUMPS link in documentation 2026-04-02 16:28:42 +02:00
sfilippone 04d7122380 Improved error handling in _bld 2026-04-02 13:53:40 +02:00
Salvatore Filippone 087dc37868 Fix handling of external packages within CBIND. 2026-03-25 16:47:30 +01:00
Salvatore Filippone 8e685aa3ae Fix includes for SuperLU and friends in configure 2026-03-24 19:05:35 +01:00
sfilippone 6dbe4c96f2 Take out amg_const.h 2026-03-24 11:44:11 +01:00
sfilippone f012a9d05e Add cpymat optional argument to hierarchy_bld 2026-03-20 13:48:53 +01:00
sfilippone 93ea03ef1c Add .gitattributes 2026-03-20 10:01:04 +01:00
sfilippone 5663188b8c Fix use of AR in configure for other platforms 2026-03-18 16:58:27 +01:00
sfilippone 8ac05fd00e Fix license 2026-03-18 16:12:44 +01:00
sfilippone c278c2a69a Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-03-18 14:51:54 +01:00
sfilippone eae53162af Fix licensing text 2026-03-18 14:46:31 +01:00
sfilippone 0b2a212523 Fix message in samples 2026-03-18 14:29:11 +01:00
sfilippone 7b255aaf6f Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-03-17 12:29:48 +01:00
sfilippone 4c9687c89b Multiple changes to CBIND and configure 2026-03-17 12:27:59 +01:00
Salvatore Filippone a1952a5bf8 Merge pull request #8 from fdrmrc/fix_C_interface
Fix missing library in C interface linking
2026-03-17 09:16:17 +01:00
Marco Feder d2070bcc05 Fix Makefile. Add psb_ext 2026-03-17 00:37:59 +01:00
sfilippone f760396403 Fix compilation for new use of MPIC 2026-03-13 14:40:00 +01:00
sfilippone 8f259367de Merge branch 'gpucinterfaces' of github.com:sfilippone/amg4psblas into gpucinterfaces 2026-01-27 14:49:48 +01:00
fdurastante 15d1386b73 Removed debug prints, fixed name of function in error strings 2026-01-27 14:46:35 +01:00
sfilippone 14d24b1d7c Inconsistent C function names 2026-01-27 12:32:46 +01:00
sfilippone 92425e2478 Merge branch 'gpucinterfaces' of github.com:sfilippone/amg4psblas into gpucinterfaces 2026-01-27 12:19:52 +01:00
sfilippone 228f4e46a9 Fix select type 2026-01-27 12:15:14 +01:00
fdurastante 492ae602f2 Added interface for smoother build. Improved options for preconditioner. There is still a memory error. 2026-01-23 14:29:21 +01:00
fdurastante 3230c70308 Standard input file for GPU experiment 2026-01-23 09:40:51 +01:00
fdurastante eb99c74fae Added L1-Jacobi smoother 2026-01-23 08:15:13 +01:00
fdurastante 1c69b2f635 Tester compiling, but misterious CUDA double free to be debugged 2026-01-22 14:35:30 +01:00
sfilippone 86731fe5bb Fix makefiles 2026-01-19 11:46:20 +01:00
fdurastante 7a5ba06622 Added interface for the allocate_wrk member function of a preconditioner. It handles also allocation on the CUDA/GPU side. 2026-01-16 14:08:36 +01:00
sfilippone d4c0428704 Add FCUDEFINES for CBIND 2026-01-16 12:24:22 +01:00
sfilippone 1a969488e3 Merge branch 'gpucinterfaces' of github.com:sfilippone/amg4psblas into gpucinterfaces 2026-01-15 15:27:53 +01:00
Salvatore Filippone 95382fbb01 Fix use of psb_stringc2f 2026-01-15 15:26:11 +01:00
sfilippone f586df1ac3 Fix use of stringc2f 2026-01-15 11:40:36 +01:00
Salvatore Filippone bab8a27962 Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-01-13 17:45:52 +01:00
sfilippone 28ccefa4bb Merge branch 'development' of github.com:sfilippone/amg4psblas into development 2026-01-13 17:43:25 +01:00
sfilippone 3f33b2ce71 Update docs 2026-01-13 17:41:47 +01:00
Salvatore Filippone 7279347055 Doc fixes 2025-12-30 19:42:01 +01:00
354 changed files with 11584 additions and 12150 deletions
+1 -1
View File
@@ -1,6 +1,6 @@
$Format:%d%n%n$
# Fall back version, probably last release:
1.2.0
1.2.1
# AMG4PSBLAS version file.
#
+8 -11
View File
@@ -13,15 +13,10 @@ endif()
# Check for the installation path for psblas
#if(NOT DEFINED PSBLAS_INSTALL_DIR)
# message(FATAL_ERROR "Please specify the path to the psblas installation directory using -DPSBLAS_INSTALL_DIR=<path>")
#endif()
message(STATUS "psblas directory is ${PSBLAS_INSTALL_DIR};;")
message(STATUS "psblas directory is ${PSBLAS_INSTALL_DIR};")
message(STATUS "PSBLAS DIRECTORY INC ${INCDIR}; MOD ${MODDIR}; LIB ${LIBDIR};")
#set(CMAKE_CXX_STANDARD 17) # Set cxx standard for the c++ part of the library
@@ -31,7 +26,7 @@ find_package(psblas REQUIRED PATHS ${PSBLAS_INSTALL_DIR})
if(NOT psblas_FOUND)
message(FATAL_ERROR "PSBLAS not found!")
else()
message(STATUS "Found PSBLAS: ${psblas_LIBRARIES}")
message(STATUS "Found PSBLAS: ${PSBLAS_LIBRARIES}")
endif()
if(CMAKE_BUILD_TYPE STREQUAL "Debug")
@@ -42,9 +37,9 @@ if(CMAKE_BUILD_TYPE STREQUAL "Debug")
message(STATUS "Fortran and CXX debug flags added: -g")
endif()
string(APPEND CMAKE_Fortran_FLAGS " -O2")
string(APPEND CMAKE_CXX_FLAGS " -O2")
message(STATUS "Fortran and CXX optimization flags added: -O2")
string(APPEND CMAKE_Fortran_FLAGS " -O2")
string(APPEND CMAKE_CXX_FLAGS " -O2")
message(STATUS "Fortran and CXX optimization flags added: -O2")
@@ -430,7 +425,9 @@ target_include_directories(amgcbind PUBLIC ${INCDIR} ${MODDIR})
target_link_libraries(amgcbind
#PUBLIC ${LAPACK_LINKER_FLAGS} ${LAPACK_LIBRARIES} ${LAPACK95_LIBRARIES}
#PUBLIC ${BLAS_LINKER_FLAGS} ${BLAS_LIBRARIES} ${BLAS95_LIBRARIES}
PUBLIC amgprec psblas::util psblas::linsolve psblas::prec psblas::ext psblas::cbind psblas::base) #TODO check actual libraries needed
PUBLIC amgprec psblas::util psblas::linsolve psblas::prec
psblas::ext psblas::cbind psblas::base)
#TODO check actual libraries needed
+1 -1
View File
@@ -7,7 +7,7 @@
(C) Copyright 2025 Salvatore Filippone
(C) Copyright 2025 Pasqua D'Ambra
(C) Copyright 2025 Fabio Durastante
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:
+1 -1
View File
@@ -71,7 +71,7 @@ EXTRALIBS=@EXTRA_LIBS@
#
AMGCDEFINES=$(SLUFLAGS) $(UMFFLAGS) $(SLUDISTFLAGS)
#$(PSBCDEFINES) $(MUMPSFLAGS)
#$(PSBCDEFINES) $(MUMPSFLAGS)
CDEFINES=$(AMGCDEFINES)
AMGFDEFINES=@AMGFDEFINES@ $(PSBFDEFINES)
FDEFINES=$(AMGFDEFINES)
+1 -1
View File
@@ -19,7 +19,7 @@ libdir:
amgobjs: mods
cd amgprec && $(MAKE) objs
cbnd: mods
cbnd: amgobjs
cd cbind && $(MAKE) objs
install: all
+4 -3
View File
@@ -58,7 +58,8 @@ MODOBJS=amg_base_prec_type.o amg_prec_type.o amg_prec_mod.o \
$(SMODOBJS) $(DMODOBJS) $(CMODOBJS) $(ZMODOBJS)
LOCAL_MODS=$(MODOBJS:.o=$(.mod))
LOCAL_MODS=$(MODOBJS:.o=$(.mod)) amg_c_l1_diag_solver$(.mod) amg_s_l1_diag_solver$(.mod) \
amg_d_l1_diag_solver$(.mod) amg_z_l1_diag_solver$(.mod)
LIBNAME=libamg_prec.a
all: mods objs impld
@@ -70,7 +71,7 @@ objs: mods impld
impld: mods
cd impl && $(MAKE)
lib: mods impld
lib: objs
cd impl && $(MAKE) lib
$(AR) $(HERE)/$(LIBNAME) $(MODOBJS)
$(RANLIB) $(HERE)/$(LIBNAME)
@@ -79,7 +80,7 @@ lib: mods impld
$(MODOBJS): $(PSBLAS_MODDIR)/$(PSBBASEMODNAME)$(.mod)
#amg_base_prec_type.o: amg_const.h
amg_base_prec_type.o: amg_config.h
amg_s_prec_type.o amg_d_prec_type.o amg_c_prec_type.o amg_z_prec_type.o : amg_base_prec_type.o
amg_prec_type.o: amg_s_prec_type.o amg_d_prec_type.o amg_c_prec_type.o amg_z_prec_type.o
amg_prec_mod.o: amg_prec_type.o amg_s_prec_mod.o amg_d_prec_mod.o amg_c_prec_mod.o amg_z_prec_mod.o
+1 -1
View File
@@ -84,7 +84,7 @@ module amg_base_prec_type
character(len=*), parameter :: amg_version_string_ = "1.2.0"
integer(psb_ipk_), parameter :: amg_version_major_ = 1
integer(psb_ipk_), parameter :: amg_version_minor_ = 2
integer(psb_ipk_), parameter :: amg_patchlevel_ = 0
integer(psb_ipk_), parameter :: amg_patchlevel_ = 1
type amg_ml_parms
integer(psb_ipk_) :: sweeps_pre, sweeps_post
+2 -2
View File
@@ -103,7 +103,7 @@ module amg_c_ainv_solver
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_ainv_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -225,7 +225,7 @@ module amg_c_ainv_solver
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_d_vect_type, psb_c_base_vect_type, psb_spk_, psb_ipk_
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
integer(psb_ipk_), intent(in) :: fillin,alg
real(psb_spk_), intent(in) :: thresh
type(psb_cspmat_type), intent(inout) :: wmat, zmat
+4 -7
View File
@@ -120,7 +120,7 @@ module amg_c_as_smoother
end interface
interface
subroutine amg_c_as_smoother_restr_v(sm,x,trans,work,info,data)
subroutine amg_c_as_smoother_restr_v(sm,x,trans,info,data)
import :: psb_cspmat_type, psb_c_vect_type, psb_c_base_vect_type, &
& psb_spk_, amg_c_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -128,7 +128,6 @@ module amg_c_as_smoother
class(amg_c_as_smoother_type), intent(inout) :: sm
type(psb_c_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_c_as_smoother_restr_v
@@ -150,7 +149,7 @@ module amg_c_as_smoother
end interface
interface
subroutine amg_c_as_smoother_prol_v(sm,x,trans,work,info,data)
subroutine amg_c_as_smoother_prol_v(sm,x,trans,info,data)
import :: psb_cspmat_type, psb_c_vect_type, psb_c_base_vect_type, &
& psb_spk_, amg_c_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -158,7 +157,6 @@ module amg_c_as_smoother
class(amg_c_as_smoother_type), intent(inout) :: sm
type(psb_c_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_c_as_smoother_prol_v
@@ -182,7 +180,7 @@ module amg_c_as_smoother
interface
subroutine amg_c_as_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_cspmat_type, psb_c_vect_type, psb_c_base_vect_type, &
& psb_spk_, amg_c_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -194,7 +192,6 @@ module amg_c_as_smoother
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -230,7 +227,7 @@ module amg_c_as_smoother
& psb_desc_type, psb_c_base_sparse_mat, psb_ipk_,&
& psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_as_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+1 -2
View File
@@ -116,7 +116,7 @@ module amg_c_base_ainv_mod
interface
subroutine amg_c_base_ainv_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_spk_,amg_c_base_ainv_solver_type, psb_c_vect_type, psb_ipk_
type(psb_desc_type), intent(in) :: desc_data
class(amg_c_base_ainv_solver_type), intent(inout) :: sv
@@ -124,7 +124,6 @@ module amg_c_base_ainv_mod
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+2 -3
View File
@@ -160,7 +160,7 @@ module amg_c_base_smoother_mod
interface
subroutine amg_c_base_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_base_smoother_type, psb_ipk_
@@ -171,7 +171,6 @@ module amg_c_base_smoother_mod
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -237,7 +236,7 @@ module amg_c_base_smoother_mod
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_base_smoother_type, psb_ipk_, psb_i_base_vect_type
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_base_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -143,7 +143,7 @@ module amg_c_base_solver_mod
interface
subroutine amg_c_base_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_base_solver_type, psb_ipk_
@@ -154,7 +154,6 @@ module amg_c_base_solver_mod
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -170,7 +169,7 @@ module amg_c_base_solver_mod
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_base_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -77,7 +77,7 @@ module amg_c_diag_solver
interface
subroutine amg_c_diag_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_diag_solver_type, psb_ipk_
@@ -87,7 +87,6 @@ module amg_c_diag_solver
type(psb_c_vect_type), intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -119,7 +118,7 @@ module amg_c_diag_solver
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -331,7 +330,7 @@ module amg_c_l1_diag_solver
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_l1_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -6
View File
@@ -106,7 +106,7 @@ module amg_c_gs_solver
interface
subroutine amg_c_gs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_c_gs_solver_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
implicit none
@@ -116,14 +116,13 @@ module amg_c_gs_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
type(psb_c_vect_type),intent(inout), optional :: initu
end subroutine amg_c_gs_solver_apply_vect
subroutine amg_c_bwgs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_c_bwgs_solver_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
implicit none
@@ -133,7 +132,6 @@ module amg_c_gs_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -181,7 +179,7 @@ module amg_c_gs_solver
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_gs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -195,7 +193,7 @@ module amg_c_gs_solver
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_bwgs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -64,7 +64,7 @@ module amg_c_id_solver
interface
subroutine amg_c_id_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_cspmat_type, psb_c_base_sparse_mat, &
& psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
& amg_c_id_solver_type, psb_ipk_
@@ -74,7 +74,6 @@ module amg_c_id_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -123,7 +122,7 @@ contains
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_id_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -101,7 +101,7 @@ module amg_c_ilu_solver
interface
subroutine amg_c_ilu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_c_ilu_solver_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_c_ilu_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -144,7 +143,7 @@ module amg_c_ilu_solver
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_ilu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -56,7 +56,7 @@ module amg_c_inner_mod
& psb_spk_, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
import :: amg_cprec_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
type(amg_cprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -80,18 +80,17 @@ module amg_c_inner_mod
complex(psb_spk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_cmlprec_aply
subroutine amg_cmlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,work,info)
subroutine amg_cmlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,info)
import :: psb_cspmat_type, psb_desc_type, &
& psb_spk_, psb_c_vect_type, psb_ipk_
import :: amg_cprec_type
implicit none
implicit none
type(psb_desc_type),intent(in) :: desc_data
type(amg_cprec_type), intent(inout) :: p
complex(psb_spk_),intent(in) :: alpha,beta
type(psb_c_vect_type),intent(inout) :: x
type(psb_c_vect_type),intent(inout) :: y
character,intent(in) :: trans
complex(psb_spk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_cmlprec_aply_vect
end interface amg_mlprec_aply
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_c_invk_solver
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_invk_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_c_invt_solver
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_invt_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -4
View File
@@ -106,7 +106,7 @@ module amg_c_jac_smoother
interface
subroutine amg_c_jac_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,trans,&
& sweeps,work,wv,info,init,initu)
& sweeps,wv,info,init,initu)
import :: psb_desc_type, amg_c_jac_smoother_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_
@@ -118,7 +118,6 @@ module amg_c_jac_smoother
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -151,7 +150,7 @@ module amg_c_jac_smoother
import :: psb_desc_type, amg_c_jac_smoother_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -274,7 +273,7 @@ module amg_c_jac_smoother
import :: psb_desc_type, amg_c_l1_jac_smoother_type, psb_c_vect_type, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_l1_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -345,6 +344,7 @@ contains
if (allocated(sm%sv)) then
call sm%sv%free(info)
if (info == psb_success_) deallocate(sm%sv,stat=info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
+3 -4
View File
@@ -101,7 +101,7 @@ module amg_c_jac_solver
interface
subroutine amg_c_jac_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_c_jac_solver_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_c_jac_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -143,7 +142,7 @@ module amg_c_jac_solver
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -160,7 +159,7 @@ module amg_c_jac_solver
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_l1_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -131,7 +131,7 @@ module amg_c_krm_solver
interface
subroutine amg_c_krm_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_c_krm_solver_type, psb_c_vect_type, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
implicit none
@@ -141,7 +141,6 @@ module amg_c_krm_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -174,7 +173,7 @@ module amg_c_krm_solver
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_krm_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -117,7 +117,7 @@ module amg_c_mumps_solver
interface
subroutine c_mumps_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_c_mumps_solver_type, psb_c_vect_type, psb_dpk_, psb_spk_, &
& psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type, psb_ipk_
implicit none
@@ -127,7 +127,6 @@ module amg_c_mumps_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -163,7 +162,7 @@ module amg_c_mumps_solver
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_mumps_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -4
View File
@@ -457,14 +457,13 @@ module amg_c_onelev_mod
integer(psb_ipk_), intent(out) :: info
complex(psb_spk_), optional :: work(:)
end subroutine amg_c_base_onelev_map_rstr_a
subroutine amg_c_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,work,vtx,vty)
subroutine amg_c_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,vtx,vty)
import
implicit none
class(amg_c_onelev_type), target, intent(inout) :: lv
complex(psb_spk_), intent(in) :: alpha, beta
type(psb_c_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
complex(psb_spk_), optional :: work(:)
type(psb_c_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_c_base_onelev_map_rstr_v
end interface
@@ -481,14 +480,13 @@ module amg_c_onelev_mod
complex(psb_spk_), optional :: work(:)
end subroutine amg_c_base_onelev_map_prol_a
subroutine amg_c_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,work,vtx,vty)
subroutine amg_c_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,vtx,vty)
import
implicit none
class(amg_c_onelev_type), target, intent(inout) :: lv
complex(psb_spk_), intent(in) :: alpha, beta
type(psb_c_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
complex(psb_spk_), optional :: work(:)
type(psb_c_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_c_base_onelev_map_prol_v
end interface
+16 -18
View File
@@ -193,7 +193,7 @@ module amg_c_prec_type
end interface
interface amg_precapply
subroutine amg_cprecaply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_cprecaply2_vect(prec,x,y,desc_data,info,trans)
import :: psb_cspmat_type, psb_desc_type, &
& psb_spk_, psb_c_vect_type, amg_cprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -202,9 +202,8 @@ module amg_c_prec_type
type(psb_c_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_spk_),intent(inout), optional, target :: work(:)
end subroutine amg_cprecaply2_vect
subroutine amg_cprecaply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_cprecaply1_vect(prec,x,desc_data,info,trans)
import :: psb_cspmat_type, psb_desc_type, &
& psb_spk_, psb_c_vect_type, amg_cprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -212,7 +211,6 @@ module amg_c_prec_type
type(psb_c_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_spk_),intent(inout), optional, target :: work(:)
end subroutine amg_cprecaply1_vect
subroutine amg_cprecaply(prec,x,y,desc_data,info,trans,work)
import :: psb_cspmat_type, psb_desc_type, psb_spk_, amg_cprec_type, psb_ipk_
@@ -311,7 +309,7 @@ module amg_c_prec_type
& psb_c_base_sparse_mat, psb_c_base_vect_type, &
& psb_i_base_vect_type, amg_cprec_type, psb_ipk_
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_cprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -323,14 +321,15 @@ module amg_c_prec_type
end interface amg_precbld
interface amg_hierarchy_bld
subroutine amg_c_hierarchy_bld(a,desc_a,prec,info)
subroutine amg_c_hierarchy_bld(a,desc_a,prec,info,cpymat)
import :: psb_cspmat_type, psb_desc_type, psb_spk_, &
& amg_cprec_type, psb_ipk_
implicit none
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_cprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
logical, intent(in), optional :: cpymat
! character, intent(in),optional :: upd
end subroutine amg_c_hierarchy_bld
end interface amg_hierarchy_bld
@@ -635,6 +634,11 @@ contains
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free(info)
if (psb_errstatus_fatal()) then
info=psb_err_internal_error_
call psb_errpush(info,name)
goto 9999
end if
end do
deallocate(prec%precv,stat=info)
end if
@@ -665,10 +669,6 @@ contains
info = psb_err_internal_error_; goto 9999
end if
!
! In the internals, do FREE on components,
! but do not deallocate them
!
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free_smoothers(info)
@@ -717,7 +717,7 @@ contains
!
! Top level methods.
!
subroutine amg_c_apply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_c_apply2_vect(prec,x,y,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_cprec_type), intent(inout) :: prec
@@ -725,7 +725,6 @@ contains
type(psb_c_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_spk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -733,7 +732,7 @@ contains
select type(prec)
type is (amg_cprec_type)
call amg_precapply(prec,x,y,desc_data,info,trans,work)
call amg_precapply(prec,x,y,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -748,14 +747,13 @@ contains
end subroutine amg_c_apply2_vect
subroutine amg_c_apply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_c_apply1_vect(prec,x,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_cprec_type), intent(inout) :: prec
type(psb_c_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_spk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -763,7 +761,7 @@ contains
select type(prec)
type is (amg_cprec_type)
call amg_precapply(prec,x,desc_data,info,trans,work)
call amg_precapply(prec,x,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -1006,7 +1004,7 @@ contains
integer(psb_ipk_), intent(out) :: info
class(psb_c_base_vect_type), intent(in), optional :: vmold
!
! In MLD the DESC optional argument is ignored, since
! In AMG the DESC optional argument is ignored, since
! the necessary info is contained in the various entries of the
! PRECV component.
type(psb_desc_type), intent(in), optional :: desc
+2 -3
View File
@@ -124,7 +124,7 @@ module amg_c_slu_solver
Implicit None
! Arguments
type(psb_cspmat_type), intent(in), target :: a
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_slu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -137,7 +137,7 @@ module amg_c_slu_solver
interface
subroutine amg_c_slu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
import amg_c_slu_solver_type
implicit none
@@ -147,7 +147,6 @@ module amg_c_slu_solver
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+493
View File
@@ -0,0 +1,493 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
!
! File: amg_c_sludist_solver_mod.f90
!
! Module: amg_c_sludist_solver_mod
!
! This module defines:
! - the amg_c_sludist_solver_type data structure containing the ingredients
! to interface with the SuperLU_Dist package.
! 1. The factorization is distributed (and thus exact)
!
!
!
module amg_c_sludist_solver
use iso_c_binding
use amg_c_base_solver_mod
#if (!defined(AMG_HAVE_SLUDIST)) || defined(PSB_IPK8)
type, extends(amg_c_base_solver_type) :: amg_c_sludist_solver_type
end type amg_c_sludist_solver_type
#else
type, extends(amg_c_base_solver_type) :: amg_c_sludist_solver_type
type(c_ptr) :: lufactors=c_null_ptr
integer(c_long_long) :: symbsize=0, numsize=0
contains
procedure, pass(sv) :: build => c_sludist_solver_bld
procedure, pass(sv) :: apply_a => c_sludist_solver_apply
procedure, pass(sv) :: apply_v => c_sludist_solver_apply_vect
procedure, pass(sv) :: free => c_sludist_solver_free
procedure, pass(sv) :: clear_data => c_sludist_solver_clear_data
procedure, pass(sv) :: descr => c_sludist_solver_descr
procedure, pass(sv) :: sizeof => c_sludist_solver_sizeof
procedure, nopass :: get_fmt => c_sludist_solver_get_fmt
procedure, nopass :: get_id => c_sludist_solver_get_id
procedure, pass(sv) :: is_global => c_sludist_solver_is_global
final :: c_sludist_solver_finalize
end type amg_c_sludist_solver_type
private :: c_sludist_solver_bld, c_sludist_solver_apply, &
& c_sludist_solver_free, c_sludist_solver_descr, &
& c_sludist_solver_sizeof, c_sludist_solver_apply_vect, &
& c_sludist_solver_get_fmt, c_sludist_solver_get_id, &
& c_sludist_solver_is_global, c_sludist_solver_clear_data
private :: c_sludist_solver_finalize
interface
function amg_csludist_fact(n,nl,nnz,ifrst, &
& values,rowptr,colind,lufactors,npr,npc) &
& bind(c,name='amg_csludist_fact') result(info)
use iso_c_binding
integer(c_int), value :: n,nl,nnz,ifrst,npr,npc
integer(c_int) :: info
integer(c_int) :: rowptr(*),colind(*)
complex(c_float_complex) :: values(*)
type(c_ptr) :: lufactors
end function amg_csludist_fact
end interface
interface
function amg_csludist_solve(itrans,n,nrhs, b, ldb, lufactors)&
& bind(c,name='amg_csludist_solve') result(info)
use iso_c_binding
integer(c_int) :: info
integer(c_int), value :: itrans,n,nrhs,ldb
complex(c_float_complex) :: b(ldb,*)
type(c_ptr), value :: lufactors
end function amg_csludist_solve
end interface
interface
function amg_csludist_free(lufactors)&
& bind(c,name='amg_csludist_free') result(info)
use iso_c_binding
integer(c_int) :: info
type(c_ptr), value :: lufactors
end function amg_csludist_free
end interface
contains
subroutine c_sludist_solver_apply(alpha,sv,x,beta,y,desc_data,&
& trans,work,info,init,initu)
use psb_base_mod
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_c_sludist_solver_type), intent(inout) :: sv
complex(psb_spk_),intent(inout) :: x(:)
complex(psb_spk_),intent(inout) :: y(:)
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
integer, intent(out) :: info
character, intent(in), optional :: init
complex(psb_spk_),intent(inout), optional :: initu(:)
integer :: n_row,n_col
complex(psb_spk_), pointer :: ww(:)
type(psb_ctxt_type) :: ctxt
integer :: np,me,i, err_act
character :: trans_
character(len=20) :: name='c_sludist_solver_apply'
call psb_erractionsave(err_act)
info = psb_success_
trans_ = psb_toupper(trans)
select case(trans_)
case('N')
case('T','C')
case default
call psb_errpush(psb_err_iarg_invalid_i_,name)
goto 9999
end select
!
! For non-iterative solvers, init and initu are ignored.
!
n_row = desc_data%get_local_rows()
n_col = desc_data%get_local_cols()
if (n_col <= size(work)) then
ww => work(1:n_col)
else
allocate(ww(n_col),stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_request_
call psb_errpush(info,name,i_err=(/n_col/),&
& a_err='complex(psb_spk_)')
goto 9999
end if
endif
if (info == psb_success_)&
& call psb_geaxpby(cone,x,czero,ww,desc_data,info)
select case(trans_)
case('N')
info = amg_csludist_solve(0,n_row,1,ww,n_row,sv%lufactors)
case('T')
info = amg_csludist_solve(1,n_row,1,ww,n_row,sv%lufactors)
case('C')
info = amg_csludist_solve(2,n_row,1,ww,n_row,sv%lufactors)
case default
call psb_errpush(psb_err_internal_error_,&
& name,a_err='Invalid TRANS in subsolve')
goto 9999
end select
if (info == psb_success_)&
& call psb_geaxpby(alpha,ww,beta,y,desc_data,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,&
& name,a_err='Error in subsolve')
goto 9999
endif
if (n_col > size(work)) then
deallocate(ww)
endif
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_sludist_solver_apply
subroutine c_sludist_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,wv,info,init,initu)
use psb_base_mod
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_c_sludist_solver_type), intent(inout) :: sv
type(psb_c_vect_type),intent(inout) :: x
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
type(psb_c_vect_type),intent(inout) :: wv(:)
integer, intent(out) :: info
character, intent(in), optional :: init
type(psb_c_vect_type),intent(inout), optional :: initu
complex(psb_spk_), target :: aux(0)
integer :: err_act
character(len=20) :: name='c_sludist_solver_apply_vect'
call psb_erractionsave(err_act)
info = psb_success_
!
! For non-iterative solvers, init and initu are ignored.
!
call x%v%sync()
call y%v%sync()
call sv%apply(alpha,x%v%v,beta,y%v%v,desc_data,trans,aux,info)
call y%v%set_host()
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_sludist_solver_apply_vect
subroutine c_sludist_solver_bld(a,desc_a,sv,info,b,amold,vmold,imold)
use psb_base_mod
Implicit None
! Arguments
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
type(psb_cspmat_type), intent(in), target, optional :: b
class(psb_c_base_sparse_mat), intent(in), optional :: amold
class(psb_c_base_vect_type), intent(in), optional :: vmold
class(psb_i_base_vect_type), intent(in), optional :: imold
! Local variables
type(psb_cspmat_type) :: atmp
type(psb_c_csr_sparse_mat) :: acsr
type(psb_ctxt_type) :: ctxt
integer(psb_lpk_), allocatable :: gia(:), gja(:)
integer(psb_lpk_) :: lfrst
integer(psb_ipk_) :: n_row,n_col, nrow_a, nztota, nglob, nzt, npr, npc
integer(psb_ipk_) :: ifrst, ibcheck
integer(psb_ipk_) :: np,me,i, err_act, debug_unit, debug_level
character(len=20) :: name='c_sludist_solver_bld', ch_err
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ctxt = desc_a%get_context()
call psb_info(ctxt, me, np)
npr = np
npc = 1
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),' start'
n_row = desc_a%get_local_rows()
n_col = desc_a%get_local_cols()
nglob = desc_a%get_global_rows()
!
! Strategy here is as follows: because a call to SLUDIST
! as a gobal solver is mostly done at the coarsest level,
! even if we start from a problem requiring 8 bytes, chances
! are that the global size will be suitable for 4 bytes
! anyway, so we hope for the best, and throw an error
! if something goes wrong.
!
if (nglob > huge(1_psb_ipk_)) then
write(0,*) me,' ',trim(name),': Error: overflow of local indices '
info=psb_err_internal_error_
call psb_errpush(info,name)
goto 9999
end if
call a%cscnv(atmp,info,type='csr')
! This in case we are dealing with AS
call psb_rwextd(n_row,atmp,info,b=b)
call atmp%mv_to(acsr)
nrow_a = acsr%get_nrows()
nztota = acsr%get_nzeros()
call psb_loc_to_glob(ione,lfrst,desc_a,info)
! Fix the entries to call C-base SuperLU
call psb_realloc(nztota,gja,info)
call psb_loc_to_glob(acsr%ja(1:nztota),gja(1:nztota), desc_a, info, iact='I')
acsr%ja(1:nztota) = gja(1:nztota)
acsr%ja(:) = acsr%ja(:) - 1
acsr%irp(:) = acsr%irp(:) - 1
ifrst = lfrst - 1
info = amg_csludist_fact(nglob,nrow_a,nztota,ifrst,&
& acsr%val,acsr%irp,acsr%ja,sv%lufactors,&
& npr,npc)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
ch_err='amg_csludist_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call acsr%free()
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),' end'
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_sludist_solver_bld
subroutine c_sludist_solver_free(sv,info)
Implicit None
! Arguments
class(amg_c_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='c_sludist_solver_free'
call psb_erractionsave(err_act)
info = 0
call sv%clear_data(info)
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_sludist_solver_free
subroutine c_sludist_solver_clear_data(sv,info)
Implicit None
! Arguments
class(amg_c_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='c_sludist_solver_clear_data'
call psb_erractionsave(err_act)
info = psb_success_
if (c_associated(sv%lufactors)) info = amg_csludist_free(sv%lufactors)
sv%lufactors = c_null_ptr
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_sludist_solver_clear_data
!
function c_sludist_solver_is_global(sv) result(val)
implicit none
class(amg_c_sludist_solver_type), intent(in) :: sv
logical :: val
val = .true.
end function c_sludist_solver_is_global
subroutine c_sludist_solver_finalize(sv)
Implicit None
! Arguments
type(amg_c_sludist_solver_type), intent(inout) :: sv
integer :: info
Integer :: err_act
character(len=20) :: name='c_sludist_solver_finalize'
call sv%free(info)
return
end subroutine c_sludist_solver_finalize
subroutine c_sludist_solver_descr(sv,info,iout,coarse,prefix)
Implicit None
! Arguments
class(amg_c_sludist_solver_type), intent(in) :: sv
integer, intent(out) :: info
integer, intent(in), optional :: iout
logical, intent(in), optional :: coarse
character(len=*), intent(in), optional :: prefix
! Local variables
integer :: err_act
type(psb_ctxt_type) :: ctxt
integer :: me, np
character(len=20), parameter :: name='amg_c_sludist_solver_descr'
integer :: iout_
character(1024) :: prefix_
call psb_erractionsave(err_act)
info = psb_success_
if (present(iout)) then
iout_ = iout
else
iout_ = psb_out_unit
endif
if (present(prefix)) then
prefix_ = prefix
else
prefix_ = ""
end if
write(iout_,*) trim(prefix_), ' SuperLU_Dist Sparse Factorization Solver. '
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_sludist_solver_descr
function c_sludist_solver_sizeof(sv) result(val)
implicit none
! Arguments
class(amg_c_sludist_solver_type), intent(in) :: sv
integer(psb_epk_) :: val
integer :: i
val = 2*psb_sizeof_ip + psb_sizeof_dp
val = val + sv%symbsize
val = val + sv%numsize
return
end function c_sludist_solver_sizeof
function c_sludist_solver_get_fmt() result(val)
implicit none
character(len=32) :: val
val = "SuperLU_Dist solver"
end function c_sludist_solver_get_fmt
function c_sludist_solver_get_id() result(val)
implicit none
integer(psb_ipk_) :: val
val = amg_sludist_
end function c_sludist_solver_get_id
#endif
end module amg_c_sludist_solver
+315
View File
@@ -0,0 +1,315 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
!
! File: amg_c_umf_solver_mod.f90
!
! Module: amg_c_umf_solver_mod
!
! This module defines:
! - the amg_c_umf_solver_type data structure containing the ingredients
! to interface with the UMFPACK package.
! 1. The factorization is restricted to the diagonal block of the
! current image.
!
module amg_c_umf_solver
use iso_c_binding
use amg_c_base_solver_mod
#if defined(PSB_IPK8)
type, extends(amg_c_base_solver_type) :: amg_c_umf_solver_type
end type amg_c_umf_solver_type
#else
type, extends(amg_c_base_solver_type) :: amg_c_umf_solver_type
type(c_ptr) :: symbolic=c_null_ptr, numeric=c_null_ptr
integer(c_long_long) :: symbsize=0, numsize=0
contains
procedure, pass(sv) :: build => amg_c_umf_solver_bld
procedure, pass(sv) :: apply_a => amg_c_umf_solver_apply
procedure, pass(sv) :: apply_v => amg_c_umf_solver_apply_vect
procedure, pass(sv) :: free => c_umf_solver_free
procedure, pass(sv) :: clear_data => c_umf_solver_clear_data
procedure, pass(sv) :: descr => c_umf_solver_descr
procedure, pass(sv) :: sizeof => c_umf_solver_sizeof
procedure, nopass :: get_fmt => c_umf_solver_get_fmt
procedure, nopass :: get_id => c_umf_solver_get_id
final :: c_umf_solver_finalize
end type amg_c_umf_solver_type
private :: c_umf_solver_free, c_umf_solver_descr, &
& c_umf_solver_sizeof, &
& c_umf_solver_get_fmt, c_umf_solver_get_id, &
& c_umf_solver_clear_data
private :: c_umf_solver_finalize
interface
function amg_cumf_fact(n,nnz,values,rowind,colptr,&
& symptr,numptr,ssize,nsize)&
& bind(c,name='amg_cumf_fact') result(info)
use iso_c_binding
integer(c_int), value :: n,nnz
integer(c_int) :: info
integer(c_long_long) :: ssize, nsize
integer(c_int) :: rowind(*),colptr(*)
complex(c_float_complex) :: values(*)
type(c_ptr) :: symptr, numptr
end function amg_cumf_fact
end interface
interface
function amg_cumf_solve(itrans,n,x, b, ldb, numptr)&
& bind(c,name='amg_cumf_solve') result(info)
use iso_c_binding
integer(c_int) :: info
integer(c_int), value :: itrans,n,ldb
complex(c_float_complex) :: x(*), b(ldb,*)
type(c_ptr), value :: numptr
end function amg_cumf_solve
end interface
interface
function amg_cumf_free(symptr, numptr)&
& bind(c,name='amg_cumf_free') result(info)
use iso_c_binding
integer(c_int) :: info
type(c_ptr), value :: symptr, numptr
end function amg_cumf_free
end interface
interface
subroutine amg_c_umf_solver_apply(alpha,sv,x,beta,y,desc_data,&
& trans,work,info,init,initu)
use psb_base_mod
import amg_c_umf_solver_type
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_c_umf_solver_type), intent(inout) :: sv
complex(psb_spk_),intent(inout) :: x(:)
complex(psb_spk_),intent(inout) :: y(:)
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_spk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
complex(psb_spk_),intent(inout), optional :: initu(:)
end subroutine amg_c_umf_solver_apply
end interface
interface
subroutine amg_c_umf_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,wv,info,init,initu)
use psb_base_mod
import amg_c_umf_solver_type
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_c_umf_solver_type), intent(inout) :: sv
type(psb_c_vect_type),intent(inout) :: x
type(psb_c_vect_type),intent(inout) :: y
complex(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
type(psb_c_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
type(psb_c_vect_type),intent(inout), optional :: initu
end subroutine amg_c_umf_solver_apply_vect
end interface
interface
subroutine amg_c_umf_solver_bld(a,desc_a,sv,info,b,amold,vmold,imold)
use psb_base_mod
import amg_c_umf_solver_type
Implicit None
! Arguments
type(psb_cspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_c_umf_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
type(psb_cspmat_type), intent(in), target, optional :: b
class(psb_c_base_sparse_mat), intent(in), optional :: amold
class(psb_c_base_vect_type), intent(in), optional :: vmold
class(psb_i_base_vect_type), intent(in), optional :: imold
end subroutine amg_c_umf_solver_bld
end interface
contains
subroutine c_umf_solver_free(sv,info)
Implicit None
! Arguments
class(amg_c_umf_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='c_umf_solver_free'
call psb_erractionsave(err_act)
call sv%clear_data(info)
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_umf_solver_free
subroutine c_umf_solver_clear_data(sv,info)
Implicit None
! Arguments
class(amg_c_umf_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='c_umf_solver_clear_data'
call psb_erractionsave(err_act)
info = 0
if (c_associated(sv%symbolic).and.c_associated(sv%numeric)) then
info = amg_cumf_free(sv%symbolic,sv%numeric)
if (info /= psb_success_) goto 9999
sv%symbolic = c_null_ptr
sv%numeric = c_null_ptr
sv%symbsize = 0
sv%numsize = 0
end if
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_umf_solver_clear_data
subroutine c_umf_solver_finalize(sv)
Implicit None
! Arguments
type(amg_c_umf_solver_type), intent(inout) :: sv
integer :: info
Integer :: err_act
character(len=20) :: name='c_umf_solver_finalize'
call sv%free(info)
return
end subroutine c_umf_solver_finalize
subroutine c_umf_solver_descr(sv,info,iout,coarse,prefix)
Implicit None
! Arguments
class(amg_c_umf_solver_type), intent(in) :: sv
integer, intent(out) :: info
integer, intent(in), optional :: iout
logical, intent(in), optional :: coarse
character(len=*), intent(in), optional :: prefix
! Local variables
integer :: err_act
character(len=20), parameter :: name='amg_c_umf_solver_descr'
integer :: iout_
character(1024) :: prefix_
call psb_erractionsave(err_act)
info = psb_success_
if (present(iout)) then
iout_ = iout
else
iout_ = psb_out_unit
endif
if (present(prefix)) then
prefix_ = prefix
else
prefix_ = ""
end if
write(iout_,*) trim(prefix_), ' UMFPACK Sparse Factorization Solver. '
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine c_umf_solver_descr
function c_umf_solver_sizeof(sv) result(val)
implicit none
! Arguments
class(amg_c_umf_solver_type), intent(in) :: sv
integer(psb_epk_) :: val
integer :: i
val = 2*psb_sizeof_lp
val = val + sv%symbsize
val = val + sv%numsize
return
end function c_umf_solver_sizeof
function c_umf_solver_get_fmt() result(val)
implicit none
character(len=32) :: val
val = "UMFPACK solver"
end function c_umf_solver_get_fmt
function c_umf_solver_get_id() result(val)
implicit none
integer(psb_ipk_) :: val
val = amg_umf_
end function c_umf_solver_get_id
#endif
end module amg_c_umf_solver
+2 -2
View File
@@ -103,7 +103,7 @@ module amg_d_ainv_solver
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_ainv_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -225,7 +225,7 @@ module amg_d_ainv_solver
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, psb_ipk_
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
integer(psb_ipk_), intent(in) :: fillin,alg
real(psb_dpk_), intent(in) :: thresh
type(psb_dspmat_type), intent(inout) :: wmat, zmat
+4 -7
View File
@@ -120,7 +120,7 @@ module amg_d_as_smoother
end interface
interface
subroutine amg_d_as_smoother_restr_v(sm,x,trans,work,info,data)
subroutine amg_d_as_smoother_restr_v(sm,x,trans,info,data)
import :: psb_dspmat_type, psb_d_vect_type, psb_d_base_vect_type, &
& psb_dpk_, amg_d_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -128,7 +128,6 @@ module amg_d_as_smoother
class(amg_d_as_smoother_type), intent(inout) :: sm
type(psb_d_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_d_as_smoother_restr_v
@@ -150,7 +149,7 @@ module amg_d_as_smoother
end interface
interface
subroutine amg_d_as_smoother_prol_v(sm,x,trans,work,info,data)
subroutine amg_d_as_smoother_prol_v(sm,x,trans,info,data)
import :: psb_dspmat_type, psb_d_vect_type, psb_d_base_vect_type, &
& psb_dpk_, amg_d_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -158,7 +157,6 @@ module amg_d_as_smoother
class(amg_d_as_smoother_type), intent(inout) :: sm
type(psb_d_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_d_as_smoother_prol_v
@@ -182,7 +180,7 @@ module amg_d_as_smoother
interface
subroutine amg_d_as_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_dspmat_type, psb_d_vect_type, psb_d_base_vect_type, &
& psb_dpk_, amg_d_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -194,7 +192,6 @@ module amg_d_as_smoother
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -230,7 +227,7 @@ module amg_d_as_smoother
& psb_desc_type, psb_d_base_sparse_mat, psb_ipk_,&
& psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_as_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+1 -2
View File
@@ -116,7 +116,7 @@ module amg_d_base_ainv_mod
interface
subroutine amg_d_base_ainv_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_dpk_,amg_d_base_ainv_solver_type, psb_d_vect_type, psb_ipk_
type(psb_desc_type), intent(in) :: desc_data
class(amg_d_base_ainv_solver_type), intent(inout) :: sv
@@ -124,7 +124,6 @@ module amg_d_base_ainv_mod
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+2 -3
View File
@@ -160,7 +160,7 @@ module amg_d_base_smoother_mod
interface
subroutine amg_d_base_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_base_smoother_type, psb_ipk_
@@ -171,7 +171,6 @@ module amg_d_base_smoother_mod
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -237,7 +236,7 @@ module amg_d_base_smoother_mod
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_base_smoother_type, psb_ipk_, psb_i_base_vect_type
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_base_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -143,7 +143,7 @@ module amg_d_base_solver_mod
interface
subroutine amg_d_base_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_base_solver_type, psb_ipk_
@@ -154,7 +154,6 @@ module amg_d_base_solver_mod
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -170,7 +169,7 @@ module amg_d_base_solver_mod
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_base_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -77,7 +77,7 @@ module amg_d_diag_solver
interface
subroutine amg_d_diag_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_diag_solver_type, psb_ipk_
@@ -87,7 +87,6 @@ module amg_d_diag_solver
type(psb_d_vect_type), intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -119,7 +118,7 @@ module amg_d_diag_solver
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -331,7 +330,7 @@ module amg_d_l1_diag_solver
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_l1_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -6
View File
@@ -106,7 +106,7 @@ module amg_d_gs_solver
interface
subroutine amg_d_gs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_d_gs_solver_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
implicit none
@@ -116,14 +116,13 @@ module amg_d_gs_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
type(psb_d_vect_type),intent(inout), optional :: initu
end subroutine amg_d_gs_solver_apply_vect
subroutine amg_d_bwgs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_d_bwgs_solver_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
implicit none
@@ -133,7 +132,6 @@ module amg_d_gs_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -181,7 +179,7 @@ module amg_d_gs_solver
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_gs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -195,7 +193,7 @@ module amg_d_gs_solver
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_bwgs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -64,7 +64,7 @@ module amg_d_id_solver
interface
subroutine amg_d_id_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_dspmat_type, psb_d_base_sparse_mat, &
& psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
& amg_d_id_solver_type, psb_ipk_
@@ -74,7 +74,6 @@ module amg_d_id_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -123,7 +122,7 @@ contains
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_id_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -101,7 +101,7 @@ module amg_d_ilu_solver
interface
subroutine amg_d_ilu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_d_ilu_solver_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_d_ilu_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -144,7 +143,7 @@ module amg_d_ilu_solver
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_ilu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -56,7 +56,7 @@ module amg_d_inner_mod
& psb_dpk_, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
import :: amg_dprec_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
type(amg_dprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -80,18 +80,17 @@ module amg_d_inner_mod
real(psb_dpk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_dmlprec_aply
subroutine amg_dmlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,work,info)
subroutine amg_dmlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,info)
import :: psb_dspmat_type, psb_desc_type, &
& psb_dpk_, psb_d_vect_type, psb_ipk_
import :: amg_dprec_type
implicit none
implicit none
type(psb_desc_type),intent(in) :: desc_data
type(amg_dprec_type), intent(inout) :: p
real(psb_dpk_),intent(in) :: alpha,beta
type(psb_d_vect_type),intent(inout) :: x
type(psb_d_vect_type),intent(inout) :: y
character,intent(in) :: trans
real(psb_dpk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_dmlprec_aply_vect
end interface amg_mlprec_aply
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_d_invk_solver
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_invk_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_d_invt_solver
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_invt_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -4
View File
@@ -106,7 +106,7 @@ module amg_d_jac_smoother
interface
subroutine amg_d_jac_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,trans,&
& sweeps,work,wv,info,init,initu)
& sweeps,wv,info,init,initu)
import :: psb_desc_type, amg_d_jac_smoother_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_
@@ -118,7 +118,6 @@ module amg_d_jac_smoother
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -151,7 +150,7 @@ module amg_d_jac_smoother
import :: psb_desc_type, amg_d_jac_smoother_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -274,7 +273,7 @@ module amg_d_jac_smoother
import :: psb_desc_type, amg_d_l1_jac_smoother_type, psb_d_vect_type, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_l1_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -345,6 +344,7 @@ contains
if (allocated(sm%sv)) then
call sm%sv%free(info)
if (info == psb_success_) deallocate(sm%sv,stat=info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
+3 -4
View File
@@ -101,7 +101,7 @@ module amg_d_jac_solver
interface
subroutine amg_d_jac_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_d_jac_solver_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_d_jac_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -143,7 +142,7 @@ module amg_d_jac_solver
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -160,7 +159,7 @@ module amg_d_jac_solver
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_l1_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -131,7 +131,7 @@ module amg_d_krm_solver
interface
subroutine amg_d_krm_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_d_krm_solver_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
implicit none
@@ -141,7 +141,6 @@ module amg_d_krm_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -174,7 +173,7 @@ module amg_d_krm_solver
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_krm_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -117,7 +117,7 @@ module amg_d_mumps_solver
interface
subroutine d_mumps_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_d_mumps_solver_type, psb_d_vect_type, psb_dpk_, psb_spk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type, psb_ipk_
implicit none
@@ -127,7 +127,6 @@ module amg_d_mumps_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -163,7 +162,7 @@ module amg_d_mumps_solver
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_mumps_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -4
View File
@@ -458,14 +458,13 @@ module amg_d_onelev_mod
integer(psb_ipk_), intent(out) :: info
real(psb_dpk_), optional :: work(:)
end subroutine amg_d_base_onelev_map_rstr_a
subroutine amg_d_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,work,vtx,vty)
subroutine amg_d_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,vtx,vty)
import
implicit none
class(amg_d_onelev_type), target, intent(inout) :: lv
real(psb_dpk_), intent(in) :: alpha, beta
type(psb_d_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
real(psb_dpk_), optional :: work(:)
type(psb_d_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_d_base_onelev_map_rstr_v
end interface
@@ -482,14 +481,13 @@ module amg_d_onelev_mod
real(psb_dpk_), optional :: work(:)
end subroutine amg_d_base_onelev_map_prol_a
subroutine amg_d_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,work,vtx,vty)
subroutine amg_d_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,vtx,vty)
import
implicit none
class(amg_d_onelev_type), target, intent(inout) :: lv
real(psb_dpk_), intent(in) :: alpha, beta
type(psb_d_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
real(psb_dpk_), optional :: work(:)
type(psb_d_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_d_base_onelev_map_prol_v
end interface
+3 -3
View File
@@ -94,7 +94,7 @@ module amg_d_poly_smoother
interface
subroutine amg_d_poly_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,trans,&
& sweeps,work,wv,info,init,initu)
& sweeps,wv,info,init,initu)
import :: psb_desc_type, amg_d_poly_smoother_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_
@@ -106,7 +106,6 @@ module amg_d_poly_smoother
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -140,7 +139,7 @@ module amg_d_poly_smoother
import :: psb_desc_type, amg_d_poly_smoother_type, psb_d_vect_type, psb_dpk_, &
& psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_poly_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -279,6 +278,7 @@ contains
if (allocated(sm%sv)) then
call sm%sv%free(info)
if (info == psb_success_) deallocate(sm%sv,stat=info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
+16 -18
View File
@@ -193,7 +193,7 @@ module amg_d_prec_type
end interface
interface amg_precapply
subroutine amg_dprecaply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_dprecaply2_vect(prec,x,y,desc_data,info,trans)
import :: psb_dspmat_type, psb_desc_type, &
& psb_dpk_, psb_d_vect_type, amg_dprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -202,9 +202,8 @@ module amg_d_prec_type
type(psb_d_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_dpk_),intent(inout), optional, target :: work(:)
end subroutine amg_dprecaply2_vect
subroutine amg_dprecaply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_dprecaply1_vect(prec,x,desc_data,info,trans)
import :: psb_dspmat_type, psb_desc_type, &
& psb_dpk_, psb_d_vect_type, amg_dprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -212,7 +211,6 @@ module amg_d_prec_type
type(psb_d_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_dpk_),intent(inout), optional, target :: work(:)
end subroutine amg_dprecaply1_vect
subroutine amg_dprecaply(prec,x,y,desc_data,info,trans,work)
import :: psb_dspmat_type, psb_desc_type, psb_dpk_, amg_dprec_type, psb_ipk_
@@ -311,7 +309,7 @@ module amg_d_prec_type
& psb_d_base_sparse_mat, psb_d_base_vect_type, &
& psb_i_base_vect_type, amg_dprec_type, psb_ipk_
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_dprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -323,14 +321,15 @@ module amg_d_prec_type
end interface amg_precbld
interface amg_hierarchy_bld
subroutine amg_d_hierarchy_bld(a,desc_a,prec,info)
subroutine amg_d_hierarchy_bld(a,desc_a,prec,info,cpymat)
import :: psb_dspmat_type, psb_desc_type, psb_dpk_, &
& amg_dprec_type, psb_ipk_
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_dprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
logical, intent(in), optional :: cpymat
! character, intent(in),optional :: upd
end subroutine amg_d_hierarchy_bld
end interface amg_hierarchy_bld
@@ -635,6 +634,11 @@ contains
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free(info)
if (psb_errstatus_fatal()) then
info=psb_err_internal_error_
call psb_errpush(info,name)
goto 9999
end if
end do
deallocate(prec%precv,stat=info)
end if
@@ -665,10 +669,6 @@ contains
info = psb_err_internal_error_; goto 9999
end if
!
! In the internals, do FREE on components,
! but do not deallocate them
!
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free_smoothers(info)
@@ -717,7 +717,7 @@ contains
!
! Top level methods.
!
subroutine amg_d_apply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_d_apply2_vect(prec,x,y,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_dprec_type), intent(inout) :: prec
@@ -725,7 +725,6 @@ contains
type(psb_d_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_dpk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -733,7 +732,7 @@ contains
select type(prec)
type is (amg_dprec_type)
call amg_precapply(prec,x,y,desc_data,info,trans,work)
call amg_precapply(prec,x,y,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -748,14 +747,13 @@ contains
end subroutine amg_d_apply2_vect
subroutine amg_d_apply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_d_apply1_vect(prec,x,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_dprec_type), intent(inout) :: prec
type(psb_d_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_dpk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -763,7 +761,7 @@ contains
select type(prec)
type is (amg_dprec_type)
call amg_precapply(prec,x,desc_data,info,trans,work)
call amg_precapply(prec,x,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -1006,7 +1004,7 @@ contains
integer(psb_ipk_), intent(out) :: info
class(psb_d_base_vect_type), intent(in), optional :: vmold
!
! In MLD the DESC optional argument is ignored, since
! In AMG the DESC optional argument is ignored, since
! the necessary info is contained in the various entries of the
! PRECV component.
type(psb_desc_type), intent(in), optional :: desc
+2 -3
View File
@@ -124,7 +124,7 @@ module amg_d_slu_solver
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_slu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -137,7 +137,7 @@ module amg_d_slu_solver
interface
subroutine amg_d_slu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
import amg_d_slu_solver_type
implicit none
@@ -147,7 +147,6 @@ module amg_d_slu_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+5 -5
View File
@@ -212,21 +212,21 @@ contains
end subroutine d_sludist_solver_apply
subroutine d_sludist_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
implicit none
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_d_sludist_solver_type), intent(inout) :: sv
type(psb_d_vect_type),intent(inout) :: x
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer, intent(out) :: info
character, intent(in), optional :: init
type(psb_d_vect_type),intent(inout), optional :: initu
real(psb_dpk_), target :: aux(0)
integer :: err_act
character(len=20) :: name='d_sludist_solver_apply_vect'
@@ -240,7 +240,7 @@ contains
call x%v%sync()
call y%v%sync()
call sv%apply(alpha,x%v%v,beta,y%v%v,desc_data,trans,work,info)
call sv%apply(alpha,x%v%v,beta,y%v%v,desc_data,trans,aux,info)
call y%v%set_host()
if (info /= 0) goto 9999
@@ -259,7 +259,7 @@ contains
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
+2 -3
View File
@@ -138,7 +138,7 @@ module amg_d_umf_solver
interface
subroutine amg_d_umf_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
import amg_d_umf_solver_type
implicit none
@@ -148,7 +148,6 @@ module amg_d_umf_solver
type(psb_d_vect_type),intent(inout) :: y
real(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_dpk_),target, intent(inout) :: work(:)
type(psb_d_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -163,7 +162,7 @@ module amg_d_umf_solver
Implicit None
! Arguments
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_d_umf_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -2
View File
@@ -103,7 +103,7 @@ module amg_s_ainv_solver
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_ainv_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -225,7 +225,7 @@ module amg_s_ainv_solver
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_d_vect_type, psb_s_base_vect_type, psb_spk_, psb_ipk_
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
integer(psb_ipk_), intent(in) :: fillin,alg
real(psb_spk_), intent(in) :: thresh
type(psb_sspmat_type), intent(inout) :: wmat, zmat
+4 -7
View File
@@ -120,7 +120,7 @@ module amg_s_as_smoother
end interface
interface
subroutine amg_s_as_smoother_restr_v(sm,x,trans,work,info,data)
subroutine amg_s_as_smoother_restr_v(sm,x,trans,info,data)
import :: psb_sspmat_type, psb_s_vect_type, psb_s_base_vect_type, &
& psb_spk_, amg_s_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -128,7 +128,6 @@ module amg_s_as_smoother
class(amg_s_as_smoother_type), intent(inout) :: sm
type(psb_s_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_s_as_smoother_restr_v
@@ -150,7 +149,7 @@ module amg_s_as_smoother
end interface
interface
subroutine amg_s_as_smoother_prol_v(sm,x,trans,work,info,data)
subroutine amg_s_as_smoother_prol_v(sm,x,trans,info,data)
import :: psb_sspmat_type, psb_s_vect_type, psb_s_base_vect_type, &
& psb_spk_, amg_s_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -158,7 +157,6 @@ module amg_s_as_smoother
class(amg_s_as_smoother_type), intent(inout) :: sm
type(psb_s_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_s_as_smoother_prol_v
@@ -182,7 +180,7 @@ module amg_s_as_smoother
interface
subroutine amg_s_as_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_sspmat_type, psb_s_vect_type, psb_s_base_vect_type, &
& psb_spk_, amg_s_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -194,7 +192,6 @@ module amg_s_as_smoother
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -230,7 +227,7 @@ module amg_s_as_smoother
& psb_desc_type, psb_s_base_sparse_mat, psb_ipk_,&
& psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_as_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+1 -2
View File
@@ -116,7 +116,7 @@ module amg_s_base_ainv_mod
interface
subroutine amg_s_base_ainv_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_spk_,amg_s_base_ainv_solver_type, psb_s_vect_type, psb_ipk_
type(psb_desc_type), intent(in) :: desc_data
class(amg_s_base_ainv_solver_type), intent(inout) :: sv
@@ -124,7 +124,6 @@ module amg_s_base_ainv_mod
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+2 -3
View File
@@ -160,7 +160,7 @@ module amg_s_base_smoother_mod
interface
subroutine amg_s_base_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_base_smoother_type, psb_ipk_
@@ -171,7 +171,6 @@ module amg_s_base_smoother_mod
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -237,7 +236,7 @@ module amg_s_base_smoother_mod
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_base_smoother_type, psb_ipk_, psb_i_base_vect_type
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_base_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -143,7 +143,7 @@ module amg_s_base_solver_mod
interface
subroutine amg_s_base_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_base_solver_type, psb_ipk_
@@ -154,7 +154,6 @@ module amg_s_base_solver_mod
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -170,7 +169,7 @@ module amg_s_base_solver_mod
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_base_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -77,7 +77,7 @@ module amg_s_diag_solver
interface
subroutine amg_s_diag_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_diag_solver_type, psb_ipk_
@@ -87,7 +87,6 @@ module amg_s_diag_solver
type(psb_s_vect_type), intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -119,7 +118,7 @@ module amg_s_diag_solver
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -331,7 +330,7 @@ module amg_s_l1_diag_solver
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_l1_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -6
View File
@@ -106,7 +106,7 @@ module amg_s_gs_solver
interface
subroutine amg_s_gs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_s_gs_solver_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
implicit none
@@ -116,14 +116,13 @@ module amg_s_gs_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
type(psb_s_vect_type),intent(inout), optional :: initu
end subroutine amg_s_gs_solver_apply_vect
subroutine amg_s_bwgs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_s_bwgs_solver_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
implicit none
@@ -133,7 +132,6 @@ module amg_s_gs_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -181,7 +179,7 @@ module amg_s_gs_solver
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_gs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -195,7 +193,7 @@ module amg_s_gs_solver
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_bwgs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -64,7 +64,7 @@ module amg_s_id_solver
interface
subroutine amg_s_id_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_sspmat_type, psb_s_base_sparse_mat, &
& psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
& amg_s_id_solver_type, psb_ipk_
@@ -74,7 +74,6 @@ module amg_s_id_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -123,7 +122,7 @@ contains
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_id_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -101,7 +101,7 @@ module amg_s_ilu_solver
interface
subroutine amg_s_ilu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_s_ilu_solver_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_s_ilu_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -144,7 +143,7 @@ module amg_s_ilu_solver
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_ilu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -56,7 +56,7 @@ module amg_s_inner_mod
& psb_spk_, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
import :: amg_sprec_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
type(amg_sprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -80,18 +80,17 @@ module amg_s_inner_mod
real(psb_spk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_smlprec_aply
subroutine amg_smlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,work,info)
subroutine amg_smlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,info)
import :: psb_sspmat_type, psb_desc_type, &
& psb_spk_, psb_s_vect_type, psb_ipk_
import :: amg_sprec_type
implicit none
implicit none
type(psb_desc_type),intent(in) :: desc_data
type(amg_sprec_type), intent(inout) :: p
real(psb_spk_),intent(in) :: alpha,beta
type(psb_s_vect_type),intent(inout) :: x
type(psb_s_vect_type),intent(inout) :: y
character,intent(in) :: trans
real(psb_spk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_smlprec_aply_vect
end interface amg_mlprec_aply
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_s_invk_solver
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_invk_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_s_invt_solver
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_invt_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -4
View File
@@ -106,7 +106,7 @@ module amg_s_jac_smoother
interface
subroutine amg_s_jac_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,trans,&
& sweeps,work,wv,info,init,initu)
& sweeps,wv,info,init,initu)
import :: psb_desc_type, amg_s_jac_smoother_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_
@@ -118,7 +118,6 @@ module amg_s_jac_smoother
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -151,7 +150,7 @@ module amg_s_jac_smoother
import :: psb_desc_type, amg_s_jac_smoother_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -274,7 +273,7 @@ module amg_s_jac_smoother
import :: psb_desc_type, amg_s_l1_jac_smoother_type, psb_s_vect_type, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_l1_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -345,6 +344,7 @@ contains
if (allocated(sm%sv)) then
call sm%sv%free(info)
if (info == psb_success_) deallocate(sm%sv,stat=info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
+3 -4
View File
@@ -101,7 +101,7 @@ module amg_s_jac_solver
interface
subroutine amg_s_jac_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_s_jac_solver_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_s_jac_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -143,7 +142,7 @@ module amg_s_jac_solver
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -160,7 +159,7 @@ module amg_s_jac_solver
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_l1_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -131,7 +131,7 @@ module amg_s_krm_solver
interface
subroutine amg_s_krm_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_s_krm_solver_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
implicit none
@@ -141,7 +141,6 @@ module amg_s_krm_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -174,7 +173,7 @@ module amg_s_krm_solver
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_krm_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -117,7 +117,7 @@ module amg_s_mumps_solver
interface
subroutine s_mumps_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_s_mumps_solver_type, psb_s_vect_type, psb_dpk_, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type, psb_ipk_
implicit none
@@ -127,7 +127,6 @@ module amg_s_mumps_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -163,7 +162,7 @@ module amg_s_mumps_solver
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_mumps_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -4
View File
@@ -458,14 +458,13 @@ module amg_s_onelev_mod
integer(psb_ipk_), intent(out) :: info
real(psb_spk_), optional :: work(:)
end subroutine amg_s_base_onelev_map_rstr_a
subroutine amg_s_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,work,vtx,vty)
subroutine amg_s_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,vtx,vty)
import
implicit none
class(amg_s_onelev_type), target, intent(inout) :: lv
real(psb_spk_), intent(in) :: alpha, beta
type(psb_s_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
real(psb_spk_), optional :: work(:)
type(psb_s_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_s_base_onelev_map_rstr_v
end interface
@@ -482,14 +481,13 @@ module amg_s_onelev_mod
real(psb_spk_), optional :: work(:)
end subroutine amg_s_base_onelev_map_prol_a
subroutine amg_s_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,work,vtx,vty)
subroutine amg_s_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,vtx,vty)
import
implicit none
class(amg_s_onelev_type), target, intent(inout) :: lv
real(psb_spk_), intent(in) :: alpha, beta
type(psb_s_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
real(psb_spk_), optional :: work(:)
type(psb_s_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_s_base_onelev_map_prol_v
end interface
+3 -3
View File
@@ -94,7 +94,7 @@ module amg_s_poly_smoother
interface
subroutine amg_s_poly_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,trans,&
& sweeps,work,wv,info,init,initu)
& sweeps,wv,info,init,initu)
import :: psb_desc_type, amg_s_poly_smoother_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_
@@ -106,7 +106,6 @@ module amg_s_poly_smoother
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -140,7 +139,7 @@ module amg_s_poly_smoother
import :: psb_desc_type, amg_s_poly_smoother_type, psb_s_vect_type, psb_spk_, &
& psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_poly_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -279,6 +278,7 @@ contains
if (allocated(sm%sv)) then
call sm%sv%free(info)
if (info == psb_success_) deallocate(sm%sv,stat=info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
+16 -18
View File
@@ -193,7 +193,7 @@ module amg_s_prec_type
end interface
interface amg_precapply
subroutine amg_sprecaply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_sprecaply2_vect(prec,x,y,desc_data,info,trans)
import :: psb_sspmat_type, psb_desc_type, &
& psb_spk_, psb_s_vect_type, amg_sprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -202,9 +202,8 @@ module amg_s_prec_type
type(psb_s_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_spk_),intent(inout), optional, target :: work(:)
end subroutine amg_sprecaply2_vect
subroutine amg_sprecaply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_sprecaply1_vect(prec,x,desc_data,info,trans)
import :: psb_sspmat_type, psb_desc_type, &
& psb_spk_, psb_s_vect_type, amg_sprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -212,7 +211,6 @@ module amg_s_prec_type
type(psb_s_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_spk_),intent(inout), optional, target :: work(:)
end subroutine amg_sprecaply1_vect
subroutine amg_sprecaply(prec,x,y,desc_data,info,trans,work)
import :: psb_sspmat_type, psb_desc_type, psb_spk_, amg_sprec_type, psb_ipk_
@@ -311,7 +309,7 @@ module amg_s_prec_type
& psb_s_base_sparse_mat, psb_s_base_vect_type, &
& psb_i_base_vect_type, amg_sprec_type, psb_ipk_
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_sprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -323,14 +321,15 @@ module amg_s_prec_type
end interface amg_precbld
interface amg_hierarchy_bld
subroutine amg_s_hierarchy_bld(a,desc_a,prec,info)
subroutine amg_s_hierarchy_bld(a,desc_a,prec,info,cpymat)
import :: psb_sspmat_type, psb_desc_type, psb_spk_, &
& amg_sprec_type, psb_ipk_
implicit none
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_sprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
logical, intent(in), optional :: cpymat
! character, intent(in),optional :: upd
end subroutine amg_s_hierarchy_bld
end interface amg_hierarchy_bld
@@ -635,6 +634,11 @@ contains
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free(info)
if (psb_errstatus_fatal()) then
info=psb_err_internal_error_
call psb_errpush(info,name)
goto 9999
end if
end do
deallocate(prec%precv,stat=info)
end if
@@ -665,10 +669,6 @@ contains
info = psb_err_internal_error_; goto 9999
end if
!
! In the internals, do FREE on components,
! but do not deallocate them
!
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free_smoothers(info)
@@ -717,7 +717,7 @@ contains
!
! Top level methods.
!
subroutine amg_s_apply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_s_apply2_vect(prec,x,y,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_sprec_type), intent(inout) :: prec
@@ -725,7 +725,6 @@ contains
type(psb_s_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_spk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -733,7 +732,7 @@ contains
select type(prec)
type is (amg_sprec_type)
call amg_precapply(prec,x,y,desc_data,info,trans,work)
call amg_precapply(prec,x,y,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -748,14 +747,13 @@ contains
end subroutine amg_s_apply2_vect
subroutine amg_s_apply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_s_apply1_vect(prec,x,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_sprec_type), intent(inout) :: prec
type(psb_s_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
real(psb_spk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -763,7 +761,7 @@ contains
select type(prec)
type is (amg_sprec_type)
call amg_precapply(prec,x,desc_data,info,trans,work)
call amg_precapply(prec,x,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -1006,7 +1004,7 @@ contains
integer(psb_ipk_), intent(out) :: info
class(psb_s_base_vect_type), intent(in), optional :: vmold
!
! In MLD the DESC optional argument is ignored, since
! In AMG the DESC optional argument is ignored, since
! the necessary info is contained in the various entries of the
! PRECV component.
type(psb_desc_type), intent(in), optional :: desc
+2 -3
View File
@@ -124,7 +124,7 @@ module amg_s_slu_solver
Implicit None
! Arguments
type(psb_sspmat_type), intent(in), target :: a
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_slu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -137,7 +137,7 @@ module amg_s_slu_solver
interface
subroutine amg_s_slu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
import amg_s_slu_solver_type
implicit none
@@ -147,7 +147,6 @@ module amg_s_slu_solver
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+493
View File
@@ -0,0 +1,493 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
!
! File: amg_s_sludist_solver_mod.f90
!
! Module: amg_s_sludist_solver_mod
!
! This module defines:
! - the amg_s_sludist_solver_type data structure containing the ingredients
! to interface with the SuperLU_Dist package.
! 1. The factorization is distributed (and thus exact)
!
!
!
module amg_s_sludist_solver
use iso_c_binding
use amg_s_base_solver_mod
#if (!defined(AMG_HAVE_SLUDIST)) || defined(PSB_IPK8)
type, extends(amg_s_base_solver_type) :: amg_s_sludist_solver_type
end type amg_s_sludist_solver_type
#else
type, extends(amg_s_base_solver_type) :: amg_s_sludist_solver_type
type(c_ptr) :: lufactors=c_null_ptr
integer(c_long_long) :: symbsize=0, numsize=0
contains
procedure, pass(sv) :: build => s_sludist_solver_bld
procedure, pass(sv) :: apply_a => s_sludist_solver_apply
procedure, pass(sv) :: apply_v => s_sludist_solver_apply_vect
procedure, pass(sv) :: free => s_sludist_solver_free
procedure, pass(sv) :: clear_data => s_sludist_solver_clear_data
procedure, pass(sv) :: descr => s_sludist_solver_descr
procedure, pass(sv) :: sizeof => s_sludist_solver_sizeof
procedure, nopass :: get_fmt => s_sludist_solver_get_fmt
procedure, nopass :: get_id => s_sludist_solver_get_id
procedure, pass(sv) :: is_global => s_sludist_solver_is_global
final :: s_sludist_solver_finalize
end type amg_s_sludist_solver_type
private :: s_sludist_solver_bld, s_sludist_solver_apply, &
& s_sludist_solver_free, s_sludist_solver_descr, &
& s_sludist_solver_sizeof, s_sludist_solver_apply_vect, &
& s_sludist_solver_get_fmt, s_sludist_solver_get_id, &
& s_sludist_solver_is_global, s_sludist_solver_clear_data
private :: s_sludist_solver_finalize
interface
function amg_ssludist_fact(n,nl,nnz,ifrst, &
& values,rowptr,colind,lufactors,npr,npc) &
& bind(c,name='amg_ssludist_fact') result(info)
use iso_c_binding
integer(c_int), value :: n,nl,nnz,ifrst,npr,npc
integer(c_int) :: info
integer(c_int) :: rowptr(*),colind(*)
real(c_float) :: values(*)
type(c_ptr) :: lufactors
end function amg_ssludist_fact
end interface
interface
function amg_ssludist_solve(itrans,n,nrhs, b, ldb, lufactors)&
& bind(c,name='amg_ssludist_solve') result(info)
use iso_c_binding
integer(c_int) :: info
integer(c_int), value :: itrans,n,nrhs,ldb
real(c_float) :: b(ldb,*)
type(c_ptr), value :: lufactors
end function amg_ssludist_solve
end interface
interface
function amg_ssludist_free(lufactors)&
& bind(c,name='amg_ssludist_free') result(info)
use iso_c_binding
integer(c_int) :: info
type(c_ptr), value :: lufactors
end function amg_ssludist_free
end interface
contains
subroutine s_sludist_solver_apply(alpha,sv,x,beta,y,desc_data,&
& trans,work,info,init,initu)
use psb_base_mod
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_s_sludist_solver_type), intent(inout) :: sv
real(psb_spk_),intent(inout) :: x(:)
real(psb_spk_),intent(inout) :: y(:)
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
integer, intent(out) :: info
character, intent(in), optional :: init
real(psb_spk_),intent(inout), optional :: initu(:)
integer :: n_row,n_col
real(psb_spk_), pointer :: ww(:)
type(psb_ctxt_type) :: ctxt
integer :: np,me,i, err_act
character :: trans_
character(len=20) :: name='s_sludist_solver_apply'
call psb_erractionsave(err_act)
info = psb_success_
trans_ = psb_toupper(trans)
select case(trans_)
case('N')
case('T','C')
case default
call psb_errpush(psb_err_iarg_invalid_i_,name)
goto 9999
end select
!
! For non-iterative solvers, init and initu are ignored.
!
n_row = desc_data%get_local_rows()
n_col = desc_data%get_local_cols()
if (n_col <= size(work)) then
ww => work(1:n_col)
else
allocate(ww(n_col),stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_request_
call psb_errpush(info,name,i_err=(/n_col/),&
& a_err='real(psb_spk_)')
goto 9999
end if
endif
if (info == psb_success_)&
& call psb_geaxpby(sone,x,szero,ww,desc_data,info)
select case(trans_)
case('N')
info = amg_ssludist_solve(0,n_row,1,ww,n_row,sv%lufactors)
case('T')
info = amg_ssludist_solve(1,n_row,1,ww,n_row,sv%lufactors)
case('C')
info = amg_ssludist_solve(2,n_row,1,ww,n_row,sv%lufactors)
case default
call psb_errpush(psb_err_internal_error_,&
& name,a_err='Invalid TRANS in subsolve')
goto 9999
end select
if (info == psb_success_)&
& call psb_geaxpby(alpha,ww,beta,y,desc_data,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,&
& name,a_err='Error in subsolve')
goto 9999
endif
if (n_col > size(work)) then
deallocate(ww)
endif
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_sludist_solver_apply
subroutine s_sludist_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,wv,info,init,initu)
use psb_base_mod
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_s_sludist_solver_type), intent(inout) :: sv
type(psb_s_vect_type),intent(inout) :: x
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
type(psb_s_vect_type),intent(inout) :: wv(:)
integer, intent(out) :: info
character, intent(in), optional :: init
type(psb_s_vect_type),intent(inout), optional :: initu
real(psb_spk_), target :: aux(0)
integer :: err_act
character(len=20) :: name='s_sludist_solver_apply_vect'
call psb_erractionsave(err_act)
info = psb_success_
!
! For non-iterative solvers, init and initu are ignored.
!
call x%v%sync()
call y%v%sync()
call sv%apply(alpha,x%v%v,beta,y%v%v,desc_data,trans,aux,info)
call y%v%set_host()
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_sludist_solver_apply_vect
subroutine s_sludist_solver_bld(a,desc_a,sv,info,b,amold,vmold,imold)
use psb_base_mod
Implicit None
! Arguments
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
type(psb_sspmat_type), intent(in), target, optional :: b
class(psb_s_base_sparse_mat), intent(in), optional :: amold
class(psb_s_base_vect_type), intent(in), optional :: vmold
class(psb_i_base_vect_type), intent(in), optional :: imold
! Local variables
type(psb_sspmat_type) :: atmp
type(psb_s_csr_sparse_mat) :: acsr
type(psb_ctxt_type) :: ctxt
integer(psb_lpk_), allocatable :: gia(:), gja(:)
integer(psb_lpk_) :: lfrst
integer(psb_ipk_) :: n_row,n_col, nrow_a, nztota, nglob, nzt, npr, npc
integer(psb_ipk_) :: ifrst, ibcheck
integer(psb_ipk_) :: np,me,i, err_act, debug_unit, debug_level
character(len=20) :: name='s_sludist_solver_bld', ch_err
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ctxt = desc_a%get_context()
call psb_info(ctxt, me, np)
npr = np
npc = 1
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),' start'
n_row = desc_a%get_local_rows()
n_col = desc_a%get_local_cols()
nglob = desc_a%get_global_rows()
!
! Strategy here is as follows: because a call to SLUDIST
! as a gobal solver is mostly done at the coarsest level,
! even if we start from a problem requiring 8 bytes, chances
! are that the global size will be suitable for 4 bytes
! anyway, so we hope for the best, and throw an error
! if something goes wrong.
!
if (nglob > huge(1_psb_ipk_)) then
write(0,*) me,' ',trim(name),': Error: overflow of local indices '
info=psb_err_internal_error_
call psb_errpush(info,name)
goto 9999
end if
call a%cscnv(atmp,info,type='csr')
! This in case we are dealing with AS
call psb_rwextd(n_row,atmp,info,b=b)
call atmp%mv_to(acsr)
nrow_a = acsr%get_nrows()
nztota = acsr%get_nzeros()
call psb_loc_to_glob(ione,lfrst,desc_a,info)
! Fix the entries to call C-base SuperLU
call psb_realloc(nztota,gja,info)
call psb_loc_to_glob(acsr%ja(1:nztota),gja(1:nztota), desc_a, info, iact='I')
acsr%ja(1:nztota) = gja(1:nztota)
acsr%ja(:) = acsr%ja(:) - 1
acsr%irp(:) = acsr%irp(:) - 1
ifrst = lfrst - 1
info = amg_ssludist_fact(nglob,nrow_a,nztota,ifrst,&
& acsr%val,acsr%irp,acsr%ja,sv%lufactors,&
& npr,npc)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
ch_err='amg_ssludist_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call acsr%free()
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),' end'
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_sludist_solver_bld
subroutine s_sludist_solver_free(sv,info)
Implicit None
! Arguments
class(amg_s_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='s_sludist_solver_free'
call psb_erractionsave(err_act)
info = 0
call sv%clear_data(info)
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_sludist_solver_free
subroutine s_sludist_solver_clear_data(sv,info)
Implicit None
! Arguments
class(amg_s_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='s_sludist_solver_clear_data'
call psb_erractionsave(err_act)
info = psb_success_
if (c_associated(sv%lufactors)) info = amg_ssludist_free(sv%lufactors)
sv%lufactors = c_null_ptr
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_sludist_solver_clear_data
!
function s_sludist_solver_is_global(sv) result(val)
implicit none
class(amg_s_sludist_solver_type), intent(in) :: sv
logical :: val
val = .true.
end function s_sludist_solver_is_global
subroutine s_sludist_solver_finalize(sv)
Implicit None
! Arguments
type(amg_s_sludist_solver_type), intent(inout) :: sv
integer :: info
Integer :: err_act
character(len=20) :: name='s_sludist_solver_finalize'
call sv%free(info)
return
end subroutine s_sludist_solver_finalize
subroutine s_sludist_solver_descr(sv,info,iout,coarse,prefix)
Implicit None
! Arguments
class(amg_s_sludist_solver_type), intent(in) :: sv
integer, intent(out) :: info
integer, intent(in), optional :: iout
logical, intent(in), optional :: coarse
character(len=*), intent(in), optional :: prefix
! Local variables
integer :: err_act
type(psb_ctxt_type) :: ctxt
integer :: me, np
character(len=20), parameter :: name='amg_s_sludist_solver_descr'
integer :: iout_
character(1024) :: prefix_
call psb_erractionsave(err_act)
info = psb_success_
if (present(iout)) then
iout_ = iout
else
iout_ = psb_out_unit
endif
if (present(prefix)) then
prefix_ = prefix
else
prefix_ = ""
end if
write(iout_,*) trim(prefix_), ' SuperLU_Dist Sparse Factorization Solver. '
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_sludist_solver_descr
function s_sludist_solver_sizeof(sv) result(val)
implicit none
! Arguments
class(amg_s_sludist_solver_type), intent(in) :: sv
integer(psb_epk_) :: val
integer :: i
val = 2*psb_sizeof_ip + psb_sizeof_dp
val = val + sv%symbsize
val = val + sv%numsize
return
end function s_sludist_solver_sizeof
function s_sludist_solver_get_fmt() result(val)
implicit none
character(len=32) :: val
val = "SuperLU_Dist solver"
end function s_sludist_solver_get_fmt
function s_sludist_solver_get_id() result(val)
implicit none
integer(psb_ipk_) :: val
val = amg_sludist_
end function s_sludist_solver_get_id
#endif
end module amg_s_sludist_solver
+315
View File
@@ -0,0 +1,315 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
!
! File: amg_s_umf_solver_mod.f90
!
! Module: amg_s_umf_solver_mod
!
! This module defines:
! - the amg_s_umf_solver_type data structure containing the ingredients
! to interface with the UMFPACK package.
! 1. The factorization is restricted to the diagonal block of the
! current image.
!
module amg_s_umf_solver
use iso_c_binding
use amg_s_base_solver_mod
#if defined(PSB_IPK8)
type, extends(amg_s_base_solver_type) :: amg_s_umf_solver_type
end type amg_s_umf_solver_type
#else
type, extends(amg_s_base_solver_type) :: amg_s_umf_solver_type
type(c_ptr) :: symbolic=c_null_ptr, numeric=c_null_ptr
integer(c_long_long) :: symbsize=0, numsize=0
contains
procedure, pass(sv) :: build => amg_s_umf_solver_bld
procedure, pass(sv) :: apply_a => amg_s_umf_solver_apply
procedure, pass(sv) :: apply_v => amg_s_umf_solver_apply_vect
procedure, pass(sv) :: free => s_umf_solver_free
procedure, pass(sv) :: clear_data => s_umf_solver_clear_data
procedure, pass(sv) :: descr => s_umf_solver_descr
procedure, pass(sv) :: sizeof => s_umf_solver_sizeof
procedure, nopass :: get_fmt => s_umf_solver_get_fmt
procedure, nopass :: get_id => s_umf_solver_get_id
final :: s_umf_solver_finalize
end type amg_s_umf_solver_type
private :: s_umf_solver_free, s_umf_solver_descr, &
& s_umf_solver_sizeof, &
& s_umf_solver_get_fmt, s_umf_solver_get_id, &
& s_umf_solver_clear_data
private :: s_umf_solver_finalize
interface
function amg_sumf_fact(n,nnz,values,rowind,colptr,&
& symptr,numptr,ssize,nsize)&
& bind(c,name='amg_sumf_fact') result(info)
use iso_c_binding
integer(c_int), value :: n,nnz
integer(c_int) :: info
integer(c_long_long) :: ssize, nsize
integer(c_int) :: rowind(*),colptr(*)
real(c_float) :: values(*)
type(c_ptr) :: symptr, numptr
end function amg_sumf_fact
end interface
interface
function amg_sumf_solve(itrans,n,x, b, ldb, numptr)&
& bind(c,name='amg_sumf_solve') result(info)
use iso_c_binding
integer(c_int) :: info
integer(c_int), value :: itrans,n,ldb
real(c_float) :: x(*), b(ldb,*)
type(c_ptr), value :: numptr
end function amg_sumf_solve
end interface
interface
function amg_sumf_free(symptr, numptr)&
& bind(c,name='amg_sumf_free') result(info)
use iso_c_binding
integer(c_int) :: info
type(c_ptr), value :: symptr, numptr
end function amg_sumf_free
end interface
interface
subroutine amg_s_umf_solver_apply(alpha,sv,x,beta,y,desc_data,&
& trans,work,info,init,initu)
use psb_base_mod
import amg_s_umf_solver_type
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_s_umf_solver_type), intent(inout) :: sv
real(psb_spk_),intent(inout) :: x(:)
real(psb_spk_),intent(inout) :: y(:)
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
real(psb_spk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
real(psb_spk_),intent(inout), optional :: initu(:)
end subroutine amg_s_umf_solver_apply
end interface
interface
subroutine amg_s_umf_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,wv,info,init,initu)
use psb_base_mod
import amg_s_umf_solver_type
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_s_umf_solver_type), intent(inout) :: sv
type(psb_s_vect_type),intent(inout) :: x
type(psb_s_vect_type),intent(inout) :: y
real(psb_spk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
type(psb_s_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
type(psb_s_vect_type),intent(inout), optional :: initu
end subroutine amg_s_umf_solver_apply_vect
end interface
interface
subroutine amg_s_umf_solver_bld(a,desc_a,sv,info,b,amold,vmold,imold)
use psb_base_mod
import amg_s_umf_solver_type
Implicit None
! Arguments
type(psb_sspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_s_umf_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
type(psb_sspmat_type), intent(in), target, optional :: b
class(psb_s_base_sparse_mat), intent(in), optional :: amold
class(psb_s_base_vect_type), intent(in), optional :: vmold
class(psb_i_base_vect_type), intent(in), optional :: imold
end subroutine amg_s_umf_solver_bld
end interface
contains
subroutine s_umf_solver_free(sv,info)
Implicit None
! Arguments
class(amg_s_umf_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='s_umf_solver_free'
call psb_erractionsave(err_act)
call sv%clear_data(info)
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_umf_solver_free
subroutine s_umf_solver_clear_data(sv,info)
Implicit None
! Arguments
class(amg_s_umf_solver_type), intent(inout) :: sv
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='s_umf_solver_clear_data'
call psb_erractionsave(err_act)
info = 0
if (c_associated(sv%symbolic).and.c_associated(sv%numeric)) then
info = amg_sumf_free(sv%symbolic,sv%numeric)
if (info /= psb_success_) goto 9999
sv%symbolic = c_null_ptr
sv%numeric = c_null_ptr
sv%symbsize = 0
sv%numsize = 0
end if
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_umf_solver_clear_data
subroutine s_umf_solver_finalize(sv)
Implicit None
! Arguments
type(amg_s_umf_solver_type), intent(inout) :: sv
integer :: info
Integer :: err_act
character(len=20) :: name='s_umf_solver_finalize'
call sv%free(info)
return
end subroutine s_umf_solver_finalize
subroutine s_umf_solver_descr(sv,info,iout,coarse,prefix)
Implicit None
! Arguments
class(amg_s_umf_solver_type), intent(in) :: sv
integer, intent(out) :: info
integer, intent(in), optional :: iout
logical, intent(in), optional :: coarse
character(len=*), intent(in), optional :: prefix
! Local variables
integer :: err_act
character(len=20), parameter :: name='amg_s_umf_solver_descr'
integer :: iout_
character(1024) :: prefix_
call psb_erractionsave(err_act)
info = psb_success_
if (present(iout)) then
iout_ = iout
else
iout_ = psb_out_unit
endif
if (present(prefix)) then
prefix_ = prefix
else
prefix_ = ""
end if
write(iout_,*) trim(prefix_), ' UMFPACK Sparse Factorization Solver. '
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine s_umf_solver_descr
function s_umf_solver_sizeof(sv) result(val)
implicit none
! Arguments
class(amg_s_umf_solver_type), intent(in) :: sv
integer(psb_epk_) :: val
integer :: i
val = 2*psb_sizeof_lp
val = val + sv%symbsize
val = val + sv%numsize
return
end function s_umf_solver_sizeof
function s_umf_solver_get_fmt() result(val)
implicit none
character(len=32) :: val
val = "UMFPACK solver"
end function s_umf_solver_get_fmt
function s_umf_solver_get_id() result(val)
implicit none
integer(psb_ipk_) :: val
val = amg_umf_
end function s_umf_solver_get_id
#endif
end module amg_s_umf_solver
+2 -2
View File
@@ -103,7 +103,7 @@ module amg_z_ainv_solver
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_ainv_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -225,7 +225,7 @@ module amg_z_ainv_solver
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_d_vect_type, psb_z_base_vect_type, psb_dpk_, psb_ipk_
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
integer(psb_ipk_), intent(in) :: fillin,alg
real(psb_dpk_), intent(in) :: thresh
type(psb_zspmat_type), intent(inout) :: wmat, zmat
+4 -7
View File
@@ -120,7 +120,7 @@ module amg_z_as_smoother
end interface
interface
subroutine amg_z_as_smoother_restr_v(sm,x,trans,work,info,data)
subroutine amg_z_as_smoother_restr_v(sm,x,trans,info,data)
import :: psb_zspmat_type, psb_z_vect_type, psb_z_base_vect_type, &
& psb_dpk_, amg_z_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -128,7 +128,6 @@ module amg_z_as_smoother
class(amg_z_as_smoother_type), intent(inout) :: sm
type(psb_z_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_z_as_smoother_restr_v
@@ -150,7 +149,7 @@ module amg_z_as_smoother
end interface
interface
subroutine amg_z_as_smoother_prol_v(sm,x,trans,work,info,data)
subroutine amg_z_as_smoother_prol_v(sm,x,trans,info,data)
import :: psb_zspmat_type, psb_z_vect_type, psb_z_base_vect_type, &
& psb_dpk_, amg_z_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -158,7 +157,6 @@ module amg_z_as_smoother
class(amg_z_as_smoother_type), intent(inout) :: sm
type(psb_z_vect_type),intent(inout) :: x
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), optional, intent(in) :: data
end subroutine amg_z_as_smoother_prol_v
@@ -182,7 +180,7 @@ module amg_z_as_smoother
interface
subroutine amg_z_as_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_zspmat_type, psb_z_vect_type, psb_z_base_vect_type, &
& psb_dpk_, amg_z_as_smoother_type, psb_epk_, &
& psb_desc_type, psb_ipk_
@@ -194,7 +192,6 @@ module amg_z_as_smoother
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -230,7 +227,7 @@ module amg_z_as_smoother
& psb_desc_type, psb_z_base_sparse_mat, psb_ipk_,&
& psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_as_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+1 -2
View File
@@ -116,7 +116,7 @@ module amg_z_base_ainv_mod
interface
subroutine amg_z_base_ainv_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_dpk_,amg_z_base_ainv_solver_type, psb_z_vect_type, psb_ipk_
type(psb_desc_type), intent(in) :: desc_data
class(amg_z_base_ainv_solver_type), intent(inout) :: sv
@@ -124,7 +124,6 @@ module amg_z_base_ainv_mod
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+2 -3
View File
@@ -160,7 +160,7 @@ module amg_z_base_smoother_mod
interface
subroutine amg_z_base_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,&
& trans,sweeps,work,wv,info,init,initu)
& trans,sweeps,wv,info,init,initu)
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_base_smoother_type, psb_ipk_
@@ -171,7 +171,6 @@ module amg_z_base_smoother_mod
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -237,7 +236,7 @@ module amg_z_base_smoother_mod
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_base_smoother_type, psb_ipk_, psb_i_base_vect_type
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_base_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -143,7 +143,7 @@ module amg_z_base_solver_mod
interface
subroutine amg_z_base_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_base_solver_type, psb_ipk_
@@ -154,7 +154,6 @@ module amg_z_base_solver_mod
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -170,7 +169,7 @@ module amg_z_base_solver_mod
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_base_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -77,7 +77,7 @@ module amg_z_diag_solver
interface
subroutine amg_z_diag_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_diag_solver_type, psb_ipk_
@@ -87,7 +87,6 @@ module amg_z_diag_solver
type(psb_z_vect_type), intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -119,7 +118,7 @@ module amg_z_diag_solver
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -331,7 +330,7 @@ module amg_z_l1_diag_solver
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_l1_diag_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -6
View File
@@ -106,7 +106,7 @@ module amg_z_gs_solver
interface
subroutine amg_z_gs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_z_gs_solver_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
implicit none
@@ -116,14 +116,13 @@ module amg_z_gs_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
type(psb_z_vect_type),intent(inout), optional :: initu
end subroutine amg_z_gs_solver_apply_vect
subroutine amg_z_bwgs_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_z_bwgs_solver_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
implicit none
@@ -133,7 +132,6 @@ module amg_z_gs_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -181,7 +179,7 @@ module amg_z_gs_solver
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_gs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -195,7 +193,7 @@ module amg_z_gs_solver
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_bwgs_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -64,7 +64,7 @@ module amg_z_id_solver
interface
subroutine amg_z_id_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, psb_zspmat_type, psb_z_base_sparse_mat, &
& psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
& amg_z_id_solver_type, psb_ipk_
@@ -74,7 +74,6 @@ module amg_z_id_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -123,7 +122,7 @@ contains
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_id_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -101,7 +101,7 @@ module amg_z_ilu_solver
interface
subroutine amg_z_ilu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_z_ilu_solver_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_z_ilu_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -144,7 +143,7 @@ module amg_z_ilu_solver
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_ilu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+3 -4
View File
@@ -56,7 +56,7 @@ module amg_z_inner_mod
& psb_dpk_, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
import :: amg_zprec_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
type(amg_zprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -80,18 +80,17 @@ module amg_z_inner_mod
complex(psb_dpk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_zmlprec_aply
subroutine amg_zmlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,work,info)
subroutine amg_zmlprec_aply_vect(alpha,p,x,beta,y,desc_data,trans,info)
import :: psb_zspmat_type, psb_desc_type, &
& psb_dpk_, psb_z_vect_type, psb_ipk_
import :: amg_zprec_type
implicit none
implicit none
type(psb_desc_type),intent(in) :: desc_data
type(amg_zprec_type), intent(inout) :: p
complex(psb_dpk_),intent(in) :: alpha,beta
type(psb_z_vect_type),intent(inout) :: x
type(psb_z_vect_type),intent(inout) :: y
character,intent(in) :: trans
complex(psb_dpk_),target :: work(:)
integer(psb_ipk_), intent(out) :: info
end subroutine amg_zmlprec_aply_vect
end interface amg_mlprec_aply
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_z_invk_solver
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_invk_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+1 -1
View File
@@ -94,7 +94,7 @@ module amg_z_invt_solver
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_invt_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+4 -4
View File
@@ -106,7 +106,7 @@ module amg_z_jac_smoother
interface
subroutine amg_z_jac_smoother_apply_vect(alpha,sm,x,beta,y,desc_data,trans,&
& sweeps,work,wv,info,init,initu)
& sweeps,wv,info,init,initu)
import :: psb_desc_type, amg_z_jac_smoother_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_
@@ -118,7 +118,6 @@ module amg_z_jac_smoother
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
integer(psb_ipk_), intent(in) :: sweeps
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -151,7 +150,7 @@ module amg_z_jac_smoother
import :: psb_desc_type, amg_z_jac_smoother_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -274,7 +273,7 @@ module amg_z_jac_smoother
import :: psb_desc_type, amg_z_l1_jac_smoother_type, psb_z_vect_type, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_l1_jac_smoother_type), intent(inout) :: sm
integer(psb_ipk_), intent(out) :: info
@@ -345,6 +344,7 @@ contains
if (allocated(sm%sv)) then
call sm%sv%free(info)
if (info == psb_success_) deallocate(sm%sv,stat=info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
+3 -4
View File
@@ -101,7 +101,7 @@ module amg_z_jac_solver
interface
subroutine amg_z_jac_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_z_jac_solver_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
implicit none
@@ -111,7 +111,6 @@ module amg_z_jac_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -143,7 +142,7 @@ module amg_z_jac_solver
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -160,7 +159,7 @@ module amg_z_jac_solver
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_l1_jac_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -131,7 +131,7 @@ module amg_z_krm_solver
interface
subroutine amg_z_krm_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_z_krm_solver_type, psb_z_vect_type, psb_dpk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
implicit none
@@ -141,7 +141,6 @@ module amg_z_krm_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -174,7 +173,7 @@ module amg_z_krm_solver
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
& psb_ipk_, psb_i_base_vect_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_krm_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -3
View File
@@ -117,7 +117,7 @@ module amg_z_mumps_solver
interface
subroutine z_mumps_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
import :: psb_desc_type, amg_z_mumps_solver_type, psb_z_vect_type, psb_dpk_, psb_spk_, &
& psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type, psb_ipk_
implicit none
@@ -127,7 +127,6 @@ module amg_z_mumps_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -163,7 +162,7 @@ module amg_z_mumps_solver
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_mumps_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
+2 -4
View File
@@ -457,14 +457,13 @@ module amg_z_onelev_mod
integer(psb_ipk_), intent(out) :: info
complex(psb_dpk_), optional :: work(:)
end subroutine amg_z_base_onelev_map_rstr_a
subroutine amg_z_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,work,vtx,vty)
subroutine amg_z_base_onelev_map_rstr_v(lv,alpha,vect_u,beta,vect_v,info,vtx,vty)
import
implicit none
class(amg_z_onelev_type), target, intent(inout) :: lv
complex(psb_dpk_), intent(in) :: alpha, beta
type(psb_z_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
complex(psb_dpk_), optional :: work(:)
type(psb_z_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_z_base_onelev_map_rstr_v
end interface
@@ -481,14 +480,13 @@ module amg_z_onelev_mod
complex(psb_dpk_), optional :: work(:)
end subroutine amg_z_base_onelev_map_prol_a
subroutine amg_z_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,work,vtx,vty)
subroutine amg_z_base_onelev_map_prol_v(lv,alpha,vect_v,beta,vect_u,info,vtx,vty)
import
implicit none
class(amg_z_onelev_type), target, intent(inout) :: lv
complex(psb_dpk_), intent(in) :: alpha, beta
type(psb_z_vect_type), intent(inout) :: vect_u, vect_v
integer(psb_ipk_), intent(out) :: info
complex(psb_dpk_), optional :: work(:)
type(psb_z_vect_type), optional, target, intent(inout) :: vtx,vty
end subroutine amg_z_base_onelev_map_prol_v
end interface
+16 -18
View File
@@ -193,7 +193,7 @@ module amg_z_prec_type
end interface
interface amg_precapply
subroutine amg_zprecaply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_zprecaply2_vect(prec,x,y,desc_data,info,trans)
import :: psb_zspmat_type, psb_desc_type, &
& psb_dpk_, psb_z_vect_type, amg_zprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -202,9 +202,8 @@ module amg_z_prec_type
type(psb_z_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_dpk_),intent(inout), optional, target :: work(:)
end subroutine amg_zprecaply2_vect
subroutine amg_zprecaply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_zprecaply1_vect(prec,x,desc_data,info,trans)
import :: psb_zspmat_type, psb_desc_type, &
& psb_dpk_, psb_z_vect_type, amg_zprec_type, psb_ipk_
type(psb_desc_type),intent(in) :: desc_data
@@ -212,7 +211,6 @@ module amg_z_prec_type
type(psb_z_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_dpk_),intent(inout), optional, target :: work(:)
end subroutine amg_zprecaply1_vect
subroutine amg_zprecaply(prec,x,y,desc_data,info,trans,work)
import :: psb_zspmat_type, psb_desc_type, psb_dpk_, amg_zprec_type, psb_ipk_
@@ -311,7 +309,7 @@ module amg_z_prec_type
& psb_z_base_sparse_mat, psb_z_base_vect_type, &
& psb_i_base_vect_type, amg_zprec_type, psb_ipk_
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_zprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
@@ -323,14 +321,15 @@ module amg_z_prec_type
end interface amg_precbld
interface amg_hierarchy_bld
subroutine amg_z_hierarchy_bld(a,desc_a,prec,info)
subroutine amg_z_hierarchy_bld(a,desc_a,prec,info,cpymat)
import :: psb_zspmat_type, psb_desc_type, psb_dpk_, &
& amg_zprec_type, psb_ipk_
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
type(psb_desc_type), intent(inout), target :: desc_a
class(amg_zprec_type), intent(inout), target :: prec
integer(psb_ipk_), intent(out) :: info
logical, intent(in), optional :: cpymat
! character, intent(in),optional :: upd
end subroutine amg_z_hierarchy_bld
end interface amg_hierarchy_bld
@@ -635,6 +634,11 @@ contains
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free(info)
if (psb_errstatus_fatal()) then
info=psb_err_internal_error_
call psb_errpush(info,name)
goto 9999
end if
end do
deallocate(prec%precv,stat=info)
end if
@@ -665,10 +669,6 @@ contains
info = psb_err_internal_error_; goto 9999
end if
!
! In the internals, do FREE on components,
! but do not deallocate them
!
if (allocated(prec%precv)) then
do i=1,size(prec%precv)
call prec%precv(i)%free_smoothers(info)
@@ -717,7 +717,7 @@ contains
!
! Top level methods.
!
subroutine amg_z_apply2_vect(prec,x,y,desc_data,info,trans,work)
subroutine amg_z_apply2_vect(prec,x,y,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_zprec_type), intent(inout) :: prec
@@ -725,7 +725,6 @@ contains
type(psb_z_vect_type),intent(inout) :: y
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_dpk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -733,7 +732,7 @@ contains
select type(prec)
type is (amg_zprec_type)
call amg_precapply(prec,x,y,desc_data,info,trans,work)
call amg_precapply(prec,x,y,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -748,14 +747,13 @@ contains
end subroutine amg_z_apply2_vect
subroutine amg_z_apply1_vect(prec,x,desc_data,info,trans,work)
subroutine amg_z_apply1_vect(prec,x,desc_data,info,trans)
implicit none
type(psb_desc_type),intent(in) :: desc_data
class(amg_zprec_type), intent(inout) :: prec
type(psb_z_vect_type),intent(inout) :: x
integer(psb_ipk_), intent(out) :: info
character(len=1), optional :: trans
complex(psb_dpk_),intent(inout), optional, target :: work(:)
Integer(psb_ipk_) :: err_act
character(len=20) :: name='d_prec_apply'
@@ -763,7 +761,7 @@ contains
select type(prec)
type is (amg_zprec_type)
call amg_precapply(prec,x,desc_data,info,trans,work)
call amg_precapply(prec,x,desc_data,info,trans)
class default
info = psb_err_missing_override_method_
call psb_errpush(info,name)
@@ -1006,7 +1004,7 @@ contains
integer(psb_ipk_), intent(out) :: info
class(psb_z_base_vect_type), intent(in), optional :: vmold
!
! In MLD the DESC optional argument is ignored, since
! In AMG the DESC optional argument is ignored, since
! the necessary info is contained in the various entries of the
! PRECV component.
type(psb_desc_type), intent(in), optional :: desc
+2 -3
View File
@@ -124,7 +124,7 @@ module amg_z_slu_solver
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_slu_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -137,7 +137,7 @@ module amg_z_slu_solver
interface
subroutine amg_z_slu_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
import amg_z_slu_solver_type
implicit none
@@ -147,7 +147,6 @@ module amg_z_slu_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
+5 -5
View File
@@ -212,21 +212,21 @@ contains
end subroutine z_sludist_solver_apply
subroutine z_sludist_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
implicit none
implicit none
type(psb_desc_type), intent(in) :: desc_data
class(amg_z_sludist_solver_type), intent(inout) :: sv
type(psb_z_vect_type),intent(inout) :: x
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer, intent(out) :: info
character, intent(in), optional :: init
type(psb_z_vect_type),intent(inout), optional :: initu
complex(psb_dpk_), target :: aux(0)
integer :: err_act
character(len=20) :: name='z_sludist_solver_apply_vect'
@@ -240,7 +240,7 @@ contains
call x%v%sync()
call y%v%sync()
call sv%apply(alpha,x%v%v,beta,y%v%v,desc_data,trans,work,info)
call sv%apply(alpha,x%v%v,beta,y%v%v,desc_data,trans,aux,info)
call y%v%set_host()
if (info /= 0) goto 9999
@@ -259,7 +259,7 @@ contains
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_sludist_solver_type), intent(inout) :: sv
integer, intent(out) :: info
+2 -3
View File
@@ -138,7 +138,7 @@ module amg_z_umf_solver
interface
subroutine amg_z_umf_solver_apply_vect(alpha,sv,x,beta,y,desc_data,&
& trans,work,wv,info,init,initu)
& trans,wv,info,init,initu)
use psb_base_mod
import amg_z_umf_solver_type
implicit none
@@ -148,7 +148,6 @@ module amg_z_umf_solver
type(psb_z_vect_type),intent(inout) :: y
complex(psb_dpk_),intent(in) :: alpha,beta
character(len=1),intent(in) :: trans
complex(psb_dpk_),target, intent(inout) :: work(:)
type(psb_z_vect_type),intent(inout) :: wv(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in), optional :: init
@@ -163,7 +162,7 @@ module amg_z_umf_solver
Implicit None
! Arguments
type(psb_zspmat_type), intent(in), target :: a
type(psb_zspmat_type), intent(inout), target :: a
Type(psb_desc_type), Intent(inout) :: desc_a
class(amg_z_umf_solver_type), intent(inout) :: sv
integer(psb_ipk_), intent(out) :: info
@@ -0,0 +1,161 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
! File: amg_c_parmatch_aggregator_mat_asb.f90
!
! Subroutine: amg_c_parmatch_aggregator_mat_asb
! Version: real
!
!
! From a given AC to final format, generating DESC_AC.
! This is quite involved, because in the context of aggregation based
! on parallel matching we are building the matrix hierarchy within BLD_TPROL
! as we go, especially if we have multiple sweeps, hence this code is called
! in two completely different contexts:
! 1. Within bld_tprol for the internal hierarchy
! 2. Outside, from amg_hierarchy_bld
! The solution we have found is for bld_tprol to copy its output
! into special components ag%ac ag%desc_ac etc so that:
! 1. if they are allocated, it means that bld_tprol has been already invoked, we are in
! amg_hierarchy_bld and we only need to copy them
! 2. If they are not allocated, we are within bld_tprol, and we need to actually
! perform the various needed steps.
!
! Arguments:
! ag - type(amg_c_parmatch_aggregator_type), input/output.
! The aggregator object
! parms - type(amg_cml_parms), input
! The aggregation parameters
! a - type(psb_cspmat_type), input.
! The sparse matrix structure containing the local part of
! the fine-level matrix.
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the fine-level matrix.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
! ilaggr - integer, dimension(:), input
! The mapping between the row indices of the coarse-level
! matrix and the row indices of the fine-level matrix.
! ilaggr(i)=j means that node i in the adjacency graph
! of the fine-level matrix is mapped onto node j in the
! adjacency graph of the coarse-level matrix. Note that the indices
! are assumed to be shifted so as to make sure the ranges on
! the various processes do not overlap.
! nlaggr - integer, dimension(:) input
! nlaggr(i) contains the aggregates held by process i.
! ac - type(psb_cspmat_type), inout
! The coarse matrix
! desc_ac - type(psb_desc_type), output.
! The communication descriptor of the fine-level matrix.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
!
! op_prol - type(psb_cspmat_type), input/output
! The tentative prolongator on input, the computed prolongator on output
!
! op_restr - type(psb_cspmat_type), input/output
! The restrictor operator; normally, it is the transpose of the prolongator.
!
! info - integer, output.
! Error code.
!
subroutine amg_c_parmatch_aggregator_inner_mat_asb(ag,parms,a,desc_a,&
& ac,desc_ac, op_prol,op_restr,info)
use psb_base_mod
use amg_base_prec_type
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_aggregator_inner_mat_asb
implicit none
class(amg_c_parmatch_aggregator_type), target, intent(inout) :: ag
type(amg_cml_parms), intent(inout) :: parms
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_cspmat_type), intent(inout) :: op_prol,op_restr
type(psb_cspmat_type), intent(inout) :: ac
type(psb_desc_type), intent(inout) :: desc_ac
integer(psb_ipk_), intent(out) :: info
!
type(psb_ctxt_type) :: ictxt
integer(psb_ipk_) :: np, me
type(psb_lc_coo_sparse_mat) :: acoo, bcoo
type(psb_lc_csr_sparse_mat) :: acsr1
integer(psb_ipk_) :: nzl, inl
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: err_act, debug_level, debug_unit
character(len=20) :: name='d_parmatch_inner_mat_asb'
character(len=80) :: aname
logical, parameter :: debug=.false., dump_prol_restr=.false.
if (psb_get_errstatus().ne.0) return
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
info = psb_success_
ictxt = desc_a%get_context()
call psb_info(ictxt,me,np)
if (debug) write(0,*) me,' ',trim(name),' Start:',&
& allocated(ag%ac),allocated(ag%desc_ac), allocated(ag%prol),allocated(ag%restr)
select case(parms%coarse_mat)
case(amg_distr_mat_)
! Do nothing, it has already been done in spmm_bld_ov.
case(amg_repl_mat_)
!
!
if (np>1) then
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='no repl coarse_mat_ here')
goto 9999
end if
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid amg_coarse_mat_')
goto 9999
end select
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine amg_c_parmatch_aggregator_inner_mat_asb
@@ -0,0 +1,203 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
! File: amg_c_parmatch_aggregator_mat_asb.f90
!
! Subroutine: amg_c_parmatch_aggregator_mat_asb
! Version: real
!
!
! From a given AC to final format, generating DESC_AC.
! This is quite involved, because in the context of aggregation based
! on parallel matching we are building the matrix hierarchy within BLD_TPROL
! as we go, especially if we have multiple sweeps, hence this code is called
! in two completely different contexts:
! 1. Within bld_tprol for the internal hierarchy
! 2. Outside, from amg_hierarchy_bld
! The solution we have found is for bld_tprol to copy its output
! into special components ag%ac ag%desc_ac etc so that:
! 1. if they are allocated, it means that bld_tprol has been already invoked, we are in
! amg_hierarchy_bld and we only need to copy them
! 2. If they are not allocated, we are within bld_tprol, and we need to actually
! perform the various needed steps.
!
! Arguments:
! ag - type(amg_c_parmatch_aggregator_type), input/output.
! The aggregator object
! parms - type(amg_cml_parms), input
! The aggregation parameters
! a - type(psb_cspmat_type), input.
! The sparse matrix structure containing the local part of
! the fine-level matrix.
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the fine-level matrix.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
! ilaggr - integer, dimension(:), input
! The mapping between the row indices of the coarse-level
! matrix and the row indices of the fine-level matrix.
! ilaggr(i)=j means that node i in the adjacency graph
! of the fine-level matrix is mapped onto node j in the
! adjacency graph of the coarse-level matrix. Note that the indices
! are assumed to be shifted so as to make sure the ranges on
! the various processes do not overlap.
! nlaggr - integer, dimension(:) input
! nlaggr(i) contains the aggregates held by process i.
! ac - type(psb_cspmat_type), inout
! The coarse matrix
! desc_ac - type(psb_desc_type), output.
! The communication descriptor of the fine-level matrix.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
!
! op_prol - type(psb_cspmat_type), input/output
! The tentative prolongator on input, the computed prolongator on output
!
! op_restr - type(psb_cspmat_type), input/output
! The restrictor operator; normally, it is the transpose of the prolongator.
!
! info - integer, output.
! Error code.
!
subroutine amg_c_parmatch_aggregator_mat_asb(ag,parms,a,desc_a,&
& ac,desc_ac, op_prol,op_restr,info)
use psb_base_mod
use amg_base_prec_type
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_aggregator_mat_asb
implicit none
class(amg_c_parmatch_aggregator_type), target, intent(inout) :: ag
type(amg_cml_parms), intent(inout) :: parms
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(inout) :: desc_a
type(psb_cspmat_type), intent(inout) :: op_prol,ac,op_restr
type(psb_desc_type), intent(inout) :: desc_ac
integer(psb_ipk_), intent(out) :: info
!
type(psb_ctxt_type) :: ctxt
integer(psb_ipk_) :: np, me
type(psb_lc_coo_sparse_mat) :: tmpcoo
type(psb_lcspmat_type) :: tmp_ac
integer(psb_ipk_) :: i_nr, i_nc, i_nl, nzl
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: err_act, debug_level, debug_unit
character(len=20) :: name='d_parmatch_mat_asb'
character(len=80) :: aname
logical, parameter :: debug=.false., dump_prol_restr=.false., dump_ac=.false.
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
info = psb_success_
ctxt = desc_a%get_context()
call psb_info(ctxt,me,np)
if (psb_get_errstatus().ne.0) then
write(0,*) me,' From:',trim(name),':',psb_get_errstatus()
return
end if
if (debug) write(0,*) me,' ',trim(name),' Start:',&
& allocated(ag%ac),allocated(ag%desc_ac), allocated(ag%prol),allocated(ag%restr)
select case(parms%coarse_mat)
case(amg_distr_mat_)
call ac%cscnv(info,type='csr')
call op_prol%cscnv(info,type='csr')
call op_restr%cscnv(info,type='csr')
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Done ac '
case(amg_repl_mat_)
!
! We are assuming here that an d matrix
! can hold all entries
!
if (desc_ac%get_global_rows() < huge(1_psb_ipk_) ) then
ntaggr = desc_ac%get_global_rows()
i_nr = ntaggr
else
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid amg_coarse_mat_')
goto 9999
end if
call op_prol%mv_to(tmpcoo)
nzl = tmpcoo%get_nzeros()
call psb_loc_to_glob(tmpcoo%ja(1:nzl),desc_ac,info,'I')
call tmpcoo%set_ncols(i_nr)
call op_prol%mv_from(tmpcoo)
call op_restr%mv_to(tmpcoo)
nzl = tmpcoo%get_nzeros()
call psb_loc_to_glob(tmpcoo%ia(1:nzl),desc_ac,info,'I')
call tmpcoo%set_nrows(i_nr)
call op_restr%mv_from(tmpcoo)
call psb_gather(tmp_ac,ac,desc_ac,info,root=-ione,&
& dupl=psb_dupl_add_,keeploc=.false.)
call tmp_ac%mv_to(tmpcoo)
call ac%mv_from(tmpcoo)
call psb_cdall(ctxt,desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(desc_ac,info)
!
! Now that we have the descriptors and the restrictor, we should
! update the W. But we don't, because REPL is only valid
! at the coarsest level, so no need to carry over.
!
if (info /= psb_success_) goto 9999
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid amg_coarse_mat_')
goto 9999
end select
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine amg_c_parmatch_aggregator_mat_asb
@@ -0,0 +1,244 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: amg_c_base_aggregator_mat_bld.f90
!
! Subroutine: amg_c_base_aggregator_mat_bld
! Version: c
!
! This routine builds the matrix associated to the current level of the
! multilevel preconditioner from the matrix associated to the previous level,
! by using the user-specified aggregation technique (therefore, it also builds the
! prolongation and restriction operators mapping the current level to the
! previous one and vice versa).
! The current level is regarded as the coarse one, while the previous as
! the fine one. This is in agreement with the fact that the routine is called,
! by amg_mlprec_bld, only on levels >=2.
! The coarse-level matrix A_C is built from a fine-level matrix A
! by using the Galerkin approach, i.e.
!
! A_C = P_C^T A P_C,
!
! where P_C is a prolongator from the coarse level to the fine one.
!
! A mapping from the nodes of the adjacency graph of A to the nodes of the
! adjacency graph of A_C has been computed by the amg_aggrmap_bld subroutine.
! The prolongator P_C is built here from this mapping, according to the
! value of p%iprcparm(amg_aggr_kind_), specified by the user through
! amg_cprecinit and amg_zprecset.
! On output from this routine the entries of AC, op_prol, op_restr
! are still in "global numbering" mode; this is fixed in the calling routine
! amg_c_lev_aggrmat_bld.
!
! Currently four different prolongators are implemented, corresponding to
! four aggregation algorithms:
! 1. un-smoothed aggregation,
! 2. smoothed aggregation,
! 3. "bizarre" aggregation.
! 4. minimum energy
! 1. The non-smoothed aggregation uses as prolongator the piecewise constant
! interpolation operator corresponding to the fine-to-coarse level mapping built
! by p%aggr%bld_tprol. This is called tentative prolongator.
! 2. The smoothed aggregation uses as prolongator the operator obtained by applying
! a damped Jacobi smoother to the tentative prolongator.
! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of AMG4PSBLAS.
! This prolongator still requires a deep analysis and testing and its use is
! not recommended.
! 4. Minimum energy aggregation
!
! For more details see
! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level
! Schwarz method, Computing, 63 (1999), 233-263.
! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based
! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007),
! 1181-1196.
! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner
! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008)
!
!
! The main structure is:
! 1. Perform sanity checks;
! 2. Compute prolongator/restrictor/AC
!
!
! Arguments:
! ag - type(amg_c_base_aggregator_type), input/output.
! The aggregator object
! parms - type(amg_cml_parms), input
! The aggregation parameters
! a - type(psb_cspmat_type), input.
! The sparse matrix structure containing the local part of
! the fine-level matrix.
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the fine-level matrix.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
! ilaggr - integer, dimension(:), input
! The mapping between the row indices of the coarse-level
! matrix and the row indices of the fine-level matrix.
! ilaggr(i)=j means that node i in the adjacency graph
! of the fine-level matrix is mapped onto node j in the
! adjacency graph of the coarse-level matrix. Note that the indices
! are assumed to be shifted so as to make sure the ranges on
! the various processes do not overlap.
! nlaggr - integer, dimension(:) input
! nlaggr(i) contains the aggregates held by process i.
! ac - type(psb_cspmat_type), output
! The coarse matrix on output
!
! op_prol - type(psb_cspmat_type), input/output
! The tentative prolongator on input, the computed prolongator on output
!
! op_restr - type(psb_cspmat_type), output
! The restrictor operator; normally, it is the transpose of the prolongator.
!
! info - integer, output.
! Error code.
!
subroutine amg_c_parmatch_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,&
& ac,desc_ac,op_prol,op_restr,t_prol,info)
use psb_base_mod
use amg_c_inner_mod
use amg_base_prec_type
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_aggregator_mat_bld
implicit none
class(amg_c_parmatch_aggregator_type), target, intent(inout) :: ag
type(amg_sml_parms), intent(inout) :: parms
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(inout) :: desc_a
integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:)
type(psb_lcspmat_type), intent(inout) :: t_prol
type(psb_cspmat_type), intent(out) :: op_prol,ac,op_restr
type(psb_desc_type), intent(inout) :: desc_ac
integer(psb_ipk_), intent(out) :: info
! Local variables
character(len=20) :: name
type(psb_ctxt_type) :: ictxt
integer(psb_ipk_) :: np, me
integer(psb_ipk_) :: err_act
integer(psb_ipk_) :: debug_level, debug_unit
type(psb_cspmat_type) :: atmp
name='c_parmatch_mat_bld'
if (psb_get_errstatus().ne.0) return
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
info = psb_success_
ictxt = desc_a%get_context()
call psb_info(ictxt,me,np)
!
! Build the coarse-level matrix from the fine-level one, starting from
! the mapping defined by amg_aggrmap_bld and applying the aggregation
! algorithm specified by
!
call clean_shortcuts(ag)
!
! When requesting smoothed aggregation we cannot use the
! unsmoothed shortcuts
!
select case (parms%aggr_prol)
case (amg_no_smooth_)
call amg_c_parmatch_unsmth_bld(parms%aggr_prol,ag,a,desc_a,&
ilaggr,nlaggr,parms,ac,desc_ac,op_prol,op_restr,&
t_prol,info)
case(amg_smooth_prol_,amg_l1_smooth_prol_)
call amg_c_parmatch_smth_bld(parms%aggr_prol,ag,a,desc_a,&
ilaggr,nlaggr,parms,ac,desc_ac,op_prol,op_restr,&
t_prol,info)
!!$ case(amg_biz_prol_)
!!$ call amg_caggrmat_biz_bld(a,desc_a,ilaggr,nlaggr, &
!!$ & parms,ac,desc_ac,op_prol,op_restr,info)
case(amg_min_energy_)
call amg_caggrmat_minnrg_bld(parms%aggr_prol,a,desc_a,&
ilaggr,nlaggr,parms,ac,desc_ac,op_prol,op_restr,&
t_prol,info)
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='Invalid aggr kind')
goto 9999
end select
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat asb')
goto 9999
end if
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
contains
subroutine clean_shortcuts(ag)
implicit none
class(amg_c_parmatch_aggregator_type), intent(inout) :: ag
integer(psb_ipk_) :: info
if (allocated(ag%prol)) then
call ag%prol%free()
deallocate(ag%prol)
end if
if (allocated(ag%restr)) then
call ag%restr%free()
deallocate(ag%restr)
end if
if (ag%unsmoothed_hierarchy) then
if (allocated(ag%ac)) call move_alloc(ag%ac, ag%rwa)
if (allocated(ag%desc_ac)) call move_alloc(ag%desc_ac,ag%rwdesc)
else
if (allocated(ag%ac)) then
call ag%ac%free()
deallocate(ag%ac)
end if
if (allocated(ag%desc_ac)) then
call ag%desc_ac%free(info)
deallocate(ag%desc_ac)
end if
end if
end subroutine clean_shortcuts
end subroutine amg_c_parmatch_aggregator_mat_bld
@@ -0,0 +1,470 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: amg_c_parmatch_aggregator_tprol.f90
!
! Subroutine: amg_c_parmatch_aggregator_tprol
! Version: real
!
!
subroutine amg_c_parmatch_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,t_prol,info)
use psb_base_mod
use amg_base_prec_type
use amg_c_inner_mod
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_aggregator_build_tprol
use iso_c_binding
implicit none
class(amg_c_parmatch_aggregator_type), target, intent(inout) :: ag
type(amg_sml_parms), intent(inout) :: parms
type(amg_saggr_data), intent(in) :: ag_data
type(psb_cspmat_type), intent(inout) :: a
type(psb_desc_type), intent(inout) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lcspmat_type), intent(out) :: t_prol
integer(psb_ipk_), intent(out) :: info
! Local variables
real(psb_spk_), allocatable :: tmpw(:), tmpwnxt(:)
integer(psb_lpk_), allocatable :: ixaggr(:), nxaggr(:), tlaggr(:), ivr(:)
type(psb_cspmat_type) :: a_tmp
integer(psb_ipk_) :: match_algorithm, n_sweeps
integer(psb_lpk_) :: target_csize
character(len=40) :: name, ch_err
character(len=80) :: fname, prefix_
type(psb_ctxt_type) :: ictxt
integer(psb_ipk_) :: np, me
integer(psb_ipk_) :: err_act, ierr
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_ipk_) :: i, j, k, nr, nc
integer(psb_lpk_) :: isz, num_pcols, nrac, ncac, lname, nz, x_sweeps, csz
integer(psb_lpk_) :: psz, sizes(4)
type(psb_c_csr_sparse_mat), target :: csr_prol, csr_pvi, csr_prod_res, acsr
type(psb_lc_csr_sparse_mat), target :: lcsr_prol
type(psb_desc_type), allocatable :: desc_acv(:)
type(psb_lc_coo_sparse_mat) :: tmpcoo, transp_coo
type(psb_cspmat_type), allocatable :: acv(:)
type(psb_cspmat_type), allocatable :: prolv(:), restrv(:)
type(psb_lcspmat_type) :: tmp_prol, tmp_pg, tmp_restr
type(psb_desc_type) :: tmp_desc_ac, tmp_desc_ax, tmp_desc_p
integer(psb_ipk_), save :: idx_mboxp=-1, idx_spmmbld=-1, idx_sweeps_mult=-1
logical, parameter :: dump=.false., do_timings=.false., debug=.false., &
& dump_prol_restr=.false.
name='c_parmatch_tprol'
ictxt = desc_a%get_context()
call psb_info(ictxt,me,np)
if (psb_get_errstatus().ne.0) then
write(0,*) me,trim(name),' Err_status :',psb_get_errstatus()
return
end if
if (debug) write(0,*) me,trim(name),' Start '
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
info = psb_success_
if ((do_timings).and.(idx_mboxp==-1)) &
& idx_mboxp = psb_get_timer_idx("PMC_TPROL: MatchBoxP")
if ((do_timings).and.(idx_spmmbld==-1)) &
& idx_spmmbld = psb_get_timer_idx("PMC_TPROL: spmm_bld")
if ((do_timings).and.(idx_sweeps_mult==-1)) &
& idx_sweeps_mult = psb_get_timer_idx("PMC_TPROL: sweeps_mult")
call amg_check_def(parms%ml_cycle,'Multilevel cycle',&
& amg_mult_ml_,is_legal_ml_cycle)
call amg_check_def(parms%par_aggr_alg,'Aggregation',&
& amg_coupled_aggr_,is_legal_coupled_par_aggr_alg)
call amg_check_def(parms%aggr_ord,'Ordering',&
& amg_aggr_ord_nat_,is_legal_ml_aggr_ord)
call amg_check_def(parms%aggr_thresh,'Aggr_Thresh',czero,is_legal_c_aggr_thrs)
match_algorithm = ag%matching_alg
n_sweeps = ag%n_sweeps
if (2**n_sweeps /= ag%orig_aggr_size) then
if (me == 0) then
write(debug_unit, *) 'Warning: AGGR_SIZE reset to value ',2**n_sweeps
end if
end if
if (ag_data%target_coarse_size > 0) then
target_csize = ag_data%target_coarse_size
else
target_csize = ag_data%min_coarse_size
end if
if (.true.) then
block
integer(psb_ipk_) :: ipv(2)
ipv(1) = target_csize
ipv(2) = n_sweeps
call psb_bcast(ictxt,ipv)
target_csize = ipv(1)
n_sweeps = ipv(2)
end block
else
call psb_bcast(ictxt,target_csize)
call psb_bcast(ictxt,n_sweeps)
end if
if (n_sweeps /= ag%n_sweeps) then
write(0,*) me,' Inconsistent N_SWEEPS ',n_sweeps,ag%n_sweeps
end if
!!$ if (me==0) write(0,*) 'Matching sweeps: ',n_sweeps
n_sweeps = max(1,n_sweeps)
if (debug) write(0,*) me,' Copies, with n_sweeps: ',n_sweeps,target_csize
if (ag%unsmoothed_hierarchy.and.allocated(ag%base_a)) then
call ag%base_a%cp_to(acsr)
if (ag%do_clean_zeros) call acsr%clean_zeros(info)
nr = acsr%get_nrows()
if (psb_size(ag%w) < nr) call ag%bld_default_w(nr)
isz = acsr%get_ncols()
call psb_realloc(isz,ixaggr,info)
if (info == psb_success_) &
& allocate(acv(0:n_sweeps), desc_acv(0:n_sweeps),&
& prolv(n_sweeps), restrv(n_sweeps),stat=info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call acv(0)%mv_from(acsr)
call ag%base_desc%clone(desc_acv(0),info)
else
call a%cp_to(acsr)
if (ag%do_clean_zeros) call acsr%clean_zeros(info)
nr = acsr%get_nrows()
if (psb_size(ag%w) < nr) call ag%bld_default_w(nr)
isz = acsr%get_ncols()
call psb_realloc(isz,ixaggr,info)
if (info == psb_success_) &
& allocate(acv(0:n_sweeps), desc_acv(0:n_sweeps),&
& prolv(n_sweeps), restrv(n_sweeps),stat=info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call acv(0)%mv_from(acsr)
call desc_a%clone(desc_acv(0),info)
end if
nrac = desc_acv(0)%get_local_rows()
ncac = desc_acv(0)%get_local_cols()
if (debug) write(0,*) me,' On input to level: ',nrac, ncac
if (allocated(ag%prol)) then
call ag%prol%free()
deallocate(ag%prol)
end if
if (allocated(ag%restr)) then
call ag%restr%free()
deallocate(ag%restr)
end if
if (dump) then
block
type(psb_lcspmat_type) :: lac
ivr = desc_acv(0)%get_global_indices(owned=.false.)
prefix_ = "input_a"
lname = len_trim(prefix_)
fname = trim(prefix_)
write(fname(lname+1:lname+9),'(a,i3.3,a)') '_p',me, '.mtx'
call acv(0)%print(fname,head='Debug aggregates')
call lac%cp_from(acv(0))
write(fname(lname+1:lname+13),'(a,i3.3,a)') '_p',me, '-glb.mtx'
call lac%print(fname,head='Debug aggregates',iv=ivr)
call lac%free()
end block
end if
call psb_geall(tmpw,desc_acv(0),info)
tmpw(1:nr) = ag%w(1:nr)
call psb_geasb(tmpw,desc_acv(0),info)
if (debug) then
call psb_barrier(ictxt)
if (me == 0) write(0,*) 'N_sweeps ',n_sweeps,nr,desc_acv(0)%is_ok(),target_csize
end if
!
! Prepare ag%ac, ag%desc_ac, ag%prol, ag%restr to enable
! shortcuts in mat_bld and mat_asb
! and ag%desc_ax which will be needed in backfix.
!
x_sweeps = -1
sweeps_loop: do i=1, n_sweeps
if (debug) then
call psb_barrier(ictxt)
if (me==0) write(0,*) me,trim(name),' Start sweeps_loop iteration:',i,' of ',n_sweeps
end if
!
! Building prol and restr because this algorithm is not decoupled
! On exit from matchbox_build_prol, prolv(i) is in global numbering
!
!
if (debug) write(0,*) me,' Into matchbox_build_prol ',info
if (do_timings) call psb_tic(idx_mboxp)
call amg_c_matchboxp_build_prol(tmpw,acv(i-1),desc_acv(i-1),ixaggr,nxaggr,tmp_prol,info,&
& symmetrize=ag%need_symmetrize,reproducible=ag%reproducible_matching)
if (do_timings) call psb_toc(idx_mboxp)
if (debug) write(0,*) me,' Out from matchbox_build_prol ',info
if (psb_errstatus_fatal()) write(0,*)me,trim(name),'Error fatal on exit bld_tprol',info
if (debug) then
call psb_barrier(ictxt)
!!$ write(0,*) name,' Call spmm_bld sweep:',i,n_sweeps
if (me==0) write(0,*) me,trim(name),' Calling spmm_bld NSW>1:',i,&
& desc_acv(i-1)%get_local_rows(),desc_acv(i-1)%get_local_cols(),&
& desc_acv(i-1)%get_global_rows()
end if
if (i == n_sweeps) call tmp_prol%clone(tmp_pg,info)
if (do_timings) call psb_tic(idx_spmmbld)
!
! On entry, prolv(i) is in global numbering,
!
call amg_c_parmatch_spmm_bld_ov(acv(i-1),desc_acv(i-1),ixaggr,nxaggr,parms,&
& acv(i),desc_acv(i), prolv(i),restrv(1),tmp_prol,info)
if (psb_errstatus_fatal()) write(0,*)me,trim(name),'Error fatal on exit from bld_ov(i)',info
if (debug) then
call psb_barrier(ictxt)
if (me==0) write(0,*) me,trim(name),' Done spmm_bld:',i
end if
if (do_timings) call psb_toc(idx_spmmbld)
! Keep a copy of prolv(i) in global numbering for the time being, will
! need it to build the final
! if (i == n_sweeps) call prolv(i)%clone(tmp_prol,info)
call ag%inner_mat_asb(parms,acv(i-1),desc_acv(i-1),&
& acv(i),desc_acv(i),prolv(i),restrv(1),info)
if (debug) then
call psb_barrier(ictxt)
if (me==0) write(0,*) me,trim(name),' Done mat_asb:',i,sum(nxaggr),target_csize,info
csz = sum(nxaggr)
call psb_bcast(ictxt,csz)
if (csz /= sum(nxaggr)) write(0,*) me,trim(name),' Mismatch matasb',&
& csz,sum(nxaggr),target_csize
end if
if (psb_errstatus_fatal()) write(0,*)me,trim(name),'Error fatal on entry to tmpwnxt 2'
!
! Fix wnxt
!
if (info == 0) call psb_geall(tmpwnxt,desc_acv(i),info)
if (info == 0) call psb_geasb(tmpwnxt,desc_acv(i),info,scratch=.true.)
if (info == 0) call psb_halo(tmpw,desc_acv(i-1),info)
!!$ write(0,*) trestr%get_nrows(),size(tmpwnxt),trestr%get_ncols(),size(tmpw)
if (info == 0) call psb_csmm(cone,restrv(1),tmpw,czero,tmpwnxt,info)
if (info /= psb_success_) then
write(0,*)me,trim(name),'Error from mat_asb/tmpw ',info
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mat_asb 2')
goto 9999
end if
if (i == 1) then
nrac = desc_acv(1)%get_local_rows()
!!$ write(0,*) 'Copying output w_nxt ',nrac
call psb_realloc(nrac,ag%w_nxt,info)
ag%w_nxt(1:nrac) = tmpwnxt(1:nrac)
!
! ILAGGR is fixed later on, but
! get a copy in case of an early exit
!
call psb_safe_ab_cpy(ixaggr,ilaggr,info)
end if
call psb_safe_ab_cpy(nxaggr,nlaggr,info)
call move_alloc(tmpwnxt,tmpw)
if (debug) then
if (csz /= sum(nlaggr)) write(0,*) me,trim(name),' Mismatch 2 matasb',&
& csz,sum(nlaggr),target_csize, info
end if
call acv(i-1)%free()
if ((sum(nlaggr) <= target_csize).or.(any(nlaggr==0))) then
x_sweeps = i
exit sweeps_loop
end if
if (debug) then
call psb_barrier(ictxt)
if (me==0) write(0,*) me,trim(name),' Done sweeps_loop iteration:',i,' of ',n_sweeps
end if
end do sweeps_loop
if (debug) then
call psb_barrier(ictxt)
if (me==0) write(0,*) me,trim(name),' Done sweeps_loop:',x_sweeps
end if
if (x_sweeps<=0) x_sweeps = n_sweeps
if (do_timings) call psb_tic(idx_sweeps_mult)
!
! Ok, now we have all the prolongators, including the last one in global numbering.
! Build the product of all prolongators. Need a tmp_desc_ax
! which is correct but most of the time overdimensioned
!
if (.not.allocated(ag%desc_ax)) allocate(ag%desc_ax)
!
block
integer(psb_ipk_) :: i, nnz
integer(psb_lpk_) :: ncol, ncsave
if (.not.allocated(ag%ac)) allocate(ag%ac)
if (.not.allocated(ag%desc_ac)) allocate(ag%desc_ac)
call desc_acv(x_sweeps)%clone(ag%desc_ac,info)
call desc_acv(x_sweeps)%free(info)
call acv(x_sweeps)%move_alloc(ag%ac,info)
if (.not.allocated(ag%prol)) allocate(ag%prol)
if (.not.allocated(ag%restr)) allocate(ag%restr)
call psb_cd_reinit(ag%desc_ac,info)
ncsave = ag%desc_ac%get_global_rows()
!
! Note: prolv(i) is already in local numbering
! because of the call to mat_asb in the loop above.
!
call prolv(x_sweeps)%mv_to(csr_prol)
if (debug) then
call psb_barrier(ictxt)
if (me == 0) write(0,*) 'Enter prolongator product loop ',x_sweeps
end if
do i=x_sweeps-1, 1, -1
call prolv(i)%mv_to(csr_pvi)
if (psb_errstatus_fatal()) write(0,*) me,' Fatal error in prolongator loop 1'
call psb_par_spspmm(csr_pvi,desc_acv(i),csr_prol,csr_prod_res,ag%desc_ac,info)
if ((info /=0).or.psb_errstatus_fatal()) write(0,*) me,' Fatal error in prolongator loop 2',info
call csr_pvi%free()
call csr_prod_res%mv_to_fmt(csr_prol,info)
if ((info /=0).or.psb_errstatus_fatal()) write(0,*) me,' Fatal error in prolongator loop 3',info
call csr_prol%set_ncols(ag%desc_ac%get_local_cols())
if ((info /=0).or.psb_errstatus_fatal()) write(0,*) me,' Fatal error in prolongator loop 4'
end do
call csr_prol%mv_to_lfmt(lcsr_prol,info)
nnz = lcsr_prol%get_nzeros()
call ag%desc_ac%l2gip(lcsr_prol%ja(1:nnz),info)
call lcsr_prol%set_ncols(ncsave)
if (debug) then
call psb_barrier(ictxt)
if (me == 0) write(0,*) 'Done prolongator product loop ',x_sweeps
end if
!
! Fix ILAGGR here by copying from CSR_PROL%JA
!
block
integer(psb_ipk_) :: nr
nr = lcsr_prol%get_nrows()
if (nnz /= nr) then
write(0,*) me,name,' Issue with prolongator? ',nr,nnz
end if
call psb_realloc(nr,ilaggr,info)
ilaggr(1:nnz) = lcsr_prol%ja(1:nnz)
end block
call tmp_prol%mv_from(lcsr_prol)
call psb_cdasb(ag%desc_ac,info)
call ag%ac%set_ncols(ag%desc_ac%get_local_cols())
end block
call tmp_prol%move_alloc(t_prol,info)
call t_prol%set_ncols(ag%desc_ac%get_local_cols())
call t_prol%set_nrows(desc_acv(0)%get_local_rows())
nrac = ag%desc_ac%get_local_rows()
ncac = ag%desc_ac%get_local_cols()
call psb_realloc(nrac,ag%w_nxt,info)
ag%w_nxt(1:nrac) = tmpw(1:nrac)
if (do_timings) call psb_toc(idx_sweeps_mult)
if (debug) then
call psb_barrier(ictxt)
if (me == 0) write(0,*) 'Out of build loop ',x_sweeps,': Output size:',sum(nlaggr)
end if
!call psb_set_debug_level(0)
if (dump) then
block
ivr = desc_acv(x_sweeps)%get_global_indices(owned=.false.)
prefix_ = "final_ac"
lname = len_trim(prefix_)
fname = trim(prefix_)
write(fname(lname+1:lname+9),'(a,i3.3,a)') '_p',me, '.mtx'
call acv(x_sweeps)%print(fname,head='Debug aggregates')
write(fname(lname+1:lname+13),'(a,i3.3,a)') '_p',me, '-glb.mtx'
call acv(x_sweeps)%print(fname,head='Debug aggregates',iv=ivr)
prefix_ = "final_tp"
lname = len_trim(prefix_)
fname = trim(prefix_)
write(fname(lname+1:lname+9),'(a,i3.3,a)') '_p',me, '.mtx'
call t_prol%print(fname,head='Tentative prolongator')
end block
end if
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='amg_bootCMatch_if')
goto 9999
end if
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine amg_c_parmatch_aggregator_build_tprol
@@ -0,0 +1,428 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: amg_caggrmat_smth_bld.F90
!
! Subroutine: amg_caggrmat_smth_bld
! Version: real
!
! This routine builds a coarse-level matrix A_C from a fine-level matrix A
! by using the Galerkin approach, i.e.
!
! A_C = P_C^T A P_C,
!
! where P_C is a prolongator from the coarse level to the fine one.
!
! The prolongator P_C is built according to a smoothed aggregation algorithm,
! i.e. it is obtained by applying a damped Jacobi smoother to the piecewise
! constant interpolation operator P corresponding to the fine-to-coarse level
! mapping built by the amg_aggrmap_bld subroutine:
!
! P_C = (I - omega*D^(-1)A) * P,
!
! where D is the diagonal matrix with main diagonal equal to the main diagonal
! of A, and omega is a suitable smoothing parameter. An estimate of the spectral
! radius of D^(-1)A, to be used in the computation of omega, is provided,
! according to the value of p%parms%aggr_omega_alg, specified by the user
! through amg_cprecinit and amg_zprecset.
!
! The coarse-level matrix A_C is distributed among the parallel processes or
! replicated on each of them, according to the value of p%parms%coarse_mat,
! specified by the user through amg_cprecinit and amg_zprecset.
! On output from this routine the entries of AC, op_prol, op_restr
! are still in "global numbering" mode; this is fixed in the calling routine
! aggregator%mat_bld.
!
!
! Arguments:
! dol1smoothing - Select between l1-Jacobi and Jacobi as smoother for the
! tentative prolongator
! a - type(psb_cspmat_type), input.
! The sparse matrix structure containing the local part of
! the fine-level matrix.
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the fine-level matrix.
! p - type(amg_c_onelev_type), input/output.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
! parms - type(amg_cml_parms), input
! Parameters controlling the choice of algorithm
! ac - type(psb_cspmat_type), output
! The coarse matrix on output
!
! ilaggr - integer, dimension(:), input
! The mapping between the row indices of the coarse-level
! matrix and the row indices of the fine-level matrix.
! ilaggr(i)=j means that node i in the adjacency graph
! of the fine-level matrix is mapped onto node j in the
! adjacency graph of the coarse-level matrix. Note that the indices
! are assumed to be shifted so as to make sure the ranges on
! the various processes do not overlap.
! nlaggr - integer, dimension(:) input
! nlaggr(i) contains the aggregates held by process i.
! op_prol - type(psb_cspmat_type), input/output
! The tentative prolongator on input, the computed prolongator on output
!
! op_restr - type(psb_cspmat_type), output
! The restrictor operator; normally, it is the transpose of the prolongator.
!
! info - integer, output.
! Error code.
!
subroutine amg_c_parmatch_smth_bld(dol1smoothing,ag,a,desc_a,ilaggr,nlaggr,&
parms,ac,desc_ac,op_prol,op_restr,t_prol,info)
use psb_base_mod
use amg_base_prec_type
use amg_c_inner_mod
use amg_c_base_aggregator_mod
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_smth_bld
implicit none
! Arguments
integer(psb_ipk_), intent(in) :: dol1smoothing
class(amg_c_parmatch_aggregator_type), target, intent(inout) :: ag
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(inout) :: desc_a
integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:)
type(amg_sml_parms), intent(inout) :: parms
type(psb_lcspmat_type), intent(inout) :: t_prol
type(psb_cspmat_type), intent(inout) :: op_prol,ac,op_restr
type(psb_desc_type), intent(inout) :: desc_ac
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, ip, &
& naggr, nzl,naggrm1,naggrp1, i, j, k, jd, icolF, nrw
integer(psb_ipk_) :: inaggr
type(psb_ctxt_type) :: ictxt
integer(psb_ipk_) :: np, me
character(len=20) :: name
type(psb_lc_coo_sparse_mat) :: tmpcoo, ac_coo, lcoo_restr
type(psb_c_coo_sparse_mat) :: coo_prol, coo_restr
type(psb_c_csr_sparse_mat) :: acsrf, csr_prol, acsr, tcsr
real(psb_spk_), allocatable :: adiag(:)
real(psb_spk_), allocatable :: arwsum(:),l1rwsum(:)
logical :: filter_mat
integer(psb_ipk_) :: debug_level, debug_unit, err_act
integer(psb_ipk_), parameter :: ncmax=16
real(psb_spk_) :: anorm, omega, tmp, dg, theta
logical, parameter :: debug_new=.false., dump_r=.false., dump_p=.false., debug=.false.
character(len=80) :: filename
logical, parameter :: do_timings=.false.
logical :: do_l1correction=.false.
integer(psb_ipk_), save :: idx_spspmm=-1, idx_phase1=-1, idx_gtrans=-1, idx_phase2=-1, idx_refine=-1, idx_phase3=-1
integer(psb_ipk_), save :: idx_cdasb=-1, idx_ptap=-1
name='amg_parmatch_smth_bld'
info=psb_success_
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
!debug_level = 2
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
nglob = desc_a%get_global_rows()
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
theta = parms%aggr_thresh
! Check if we have to perform l1-Jacobi or Jacobi as smoother
if(dol1smoothing.eq.amg_l1_smooth_prol_) do_l1correction=.true.
!write(0,*) me,' ',trim(name),' Start ',idx_spspmm
if ((do_timings).and.(idx_spspmm==-1)) &
& idx_spspmm = psb_get_timer_idx("PMC_SMTH_BLD: par_spspmm")
if ((do_timings).and.(idx_phase1==-1)) &
& idx_phase1 = psb_get_timer_idx("PMC_SMTH_BLD: phase1 ")
if ((do_timings).and.(idx_phase2==-1)) &
& idx_phase2 = psb_get_timer_idx("PMC_SMTH_BLD: phase2 ")
if ((do_timings).and.(idx_phase3==-1)) &
& idx_phase3 = psb_get_timer_idx("PMC_SMTH_BLD: phase3 ")
if ((do_timings).and.(idx_gtrans==-1)) &
& idx_gtrans = psb_get_timer_idx("PMC_SMTH_BLD: gtrans ")
if ((do_timings).and.(idx_refine==-1)) &
& idx_refine = psb_get_timer_idx("PMC_SMTH_BLD: refine ")
if ((do_timings).and.(idx_cdasb==-1)) &
& idx_cdasb = psb_get_timer_idx("PMC_SMTH_BLD: cdasb ")
if ((do_timings).and.(idx_ptap==-1)) &
& idx_ptap = psb_get_timer_idx("PMC_SMTH_BLD: ptap_bld ")
if (do_timings) call psb_tic(idx_phase1)
naggr = nlaggr(me+1)
ntaggr = sum(nlaggr)
naggrm1 = sum(nlaggr(1:me))
naggrp1 = sum(nlaggr(1:me+1))
filter_mat = (parms%aggr_filter == amg_filter_mat_)
!
! naggr: number of local aggregates
! nrow: local rows.
!
if (dump_p) then
block
integer(psb_lpk_), allocatable :: ivr(:), ivc(:)
integer(psb_lpk_) :: i
character(len=132) :: aname
write(0,*) me,' ',trim(name),' Dumping inp_prol/restr'
write(aname,'(a,i0,a,i0,a)') 'tprol-',desc_a%get_global_rows(),'-p',me,'.mtx'
call t_prol%print(fname=aname,head='Test ')
end block
end if
if (do_timings) call psb_tic(idx_refine)
! Get the diagonal D
adiag = a%get_diag(info)
if (info == psb_success_) &
& call psb_realloc(ncol,adiag,info)
if (info == psb_success_) &
& call psb_halo(adiag,desc_a,info)
if (info == psb_success_) call a%cp_to(acsr)
! Get the l1-diagonal of D
if (do_l1correction) then
allocate(l1rwsum(nrow))
call acsr%arwsum(l1rwsum)
if (info == psb_success_) &
& call psb_realloc(ncol,l1rwsum,info)
if (info == psb_success_) &
& call psb_halo(l1rwsum,desc_a,info)
! \tilde{D}_{i,i} = \sum_{j \ne i} |a_{i,j}|
do i=1,size(adiag)
adiag(i) = adiag(i) + l1rwsum(i) - abs(adiag(i))
end do
end if
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='sp_getdiag')
goto 9999
end if
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& ' Initial copies done.'
call acsr%cp_to_fmt(acsrf,info)
if (filter_mat) then
!
! Build the filtered matrix Af from A
!
do i=1, nrow
tmp = dzero
jd = -1
do j=acsrf%irp(i),acsrf%irp(i+1)-1
if (acsrf%ja(j) == i) jd = j
if (abs(acsrf%val(j)) < theta*sqrt(abs(adiag(i)*adiag(acsrf%ja(j))))) then
tmp=tmp+acsrf%val(j)
acsrf%val(j)=dzero
endif
enddo
if (jd == -1) then
write(0,*) name,': Warning: there is no diagonal element', i
else
acsrf%val(jd)=acsrf%val(jd)-tmp
end if
enddo
! Take out zeroed terms
call acsrf%clean_zeros(info)
end if
do i=1,size(adiag)
if (adiag(i) /= dzero) then
adiag(i) = done / adiag(i)
else
adiag(i) = done
end if
end do
if (do_timings) call psb_toc(idx_refine)
if (parms%aggr_omega_alg == amg_eig_est_) then
if (do_l1correction) then
! For l1-Jacobi this can be estimated with 1
parms%aggr_omega_val = done
else if (parms%aggr_eig == amg_max_norm_) then
allocate(arwsum(nrow))
call acsr%arwsum(arwsum)
anorm = maxval(abs(adiag(1:nrow)*arwsum(1:nrow)))
call psb_amx(ictxt,anorm)
omega = 4.d0/(3.d0*anorm)
parms%aggr_omega_val = omega
else
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid amg_aggr_eig_')
goto 9999
end if
else if (parms%aggr_omega_alg == amg_user_choice_) then
omega = parms%aggr_omega_val
else if (parms%aggr_omega_alg /= amg_user_choice_) then
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid amg_aggr_omega_alg_')
goto 9999
end if
call acsrf%scal(adiag,info)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& ' Filtering and scaling done.',info
if (info /= psb_success_) goto 9999
inaggr = naggr
call t_prol%cp_to(tmpcoo)
call psb_cdall(ictxt,desc_ac,info,nl=inaggr)
nzl = tmpcoo%get_nzeros()
call desc_ac%indxmap%g2lip_ins(tmpcoo%ja(1:nzl),info)
call tmpcoo%set_ncols(desc_ac%get_local_cols())
call tmpcoo%mv_to_ifmt(tcsr,info)
!
! Build the smoothed prolongator using either A or Af
! csr_prol = (I-w*D*A) Prol csr_prol = (I-w*D*Af) Prol
! This is always done through the variable acsrf which
! is a bit less readable, but saves space and one extra matrix copy
!
call omega_smooth(omega,acsrf)
if (do_timings) call psb_toc(idx_phase1)
if (do_timings) call psb_tic(idx_spspmm)
call psb_par_spspmm(acsrf,desc_a,tcsr,csr_prol,desc_ac,info)
call tcsr%free()
if (do_timings) call psb_toc(idx_spspmm)
if (do_timings) call psb_tic(idx_phase2)
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='spspmm 1')
goto 9999
end if
!
! Now that we have the smoothed prolongator, we can
! compute the triple product.
!
if (do_timings) call psb_tic(idx_cdasb)
call psb_cdasb(desc_ac,info)
if (do_timings) call psb_toc(idx_cdasb)
call psb_cd_reinit(desc_ac,info)
call csr_prol%mv_to_coo(coo_prol,info)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Done SPSPMM 1'
if (do_timings) call psb_tic(idx_ptap)
if (.not.allocated(ag%desc_ax)) allocate(ag%desc_ax)
call amg_ptap_bld(acsr,desc_a,nlaggr,parms,ac,&
& coo_prol,desc_ac,coo_restr,info,desc_ax=ag%desc_ax)
if (do_timings) call psb_toc(idx_ptap)
call op_prol%mv_from(coo_prol)
call op_restr%mv_from(coo_restr)
if (debug) write(0,*) me,' ',trim(name),' After mv_from',psb_get_errstatus()
if (debug) write(0,*) me,' ',trim(name),' ',ac%get_fmt(),ac%get_nrows(),ac%get_ncols(),ac%get_nzeros(),naggr,ntaggr
if (dump_r) then
block
integer(psb_lpk_), allocatable :: ivr(:), ivc(:)
integer(psb_lpk_) :: i
character(len=132) :: aname
type(psb_lcspmat_type) :: aglob
type(psb_cspmat_type) :: atmp
write(0,*) me,' ',trim(name),' Dumping prol/restr'
ivc = [(i,i=1,desc_a%get_local_cols())]
call desc_a%l2gip(ivc,info)
ivr = [(i,i=1,desc_ac%get_local_cols())]
call desc_ac%l2gip(ivr,info)
write(aname,'(a,i0,a,i0,a)') 'restr-',desc_ac%get_global_rows(),'-p',me,'.mtx'
call op_restr%print(fname=aname,head='Test ',ivc=ivc)
end block
end if
if (allocated(l1rwsum)) deallocate(l1rwsum)
if (do_timings) call psb_toc(idx_phase2)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Done smooth_aggregate '
call psb_erractionrestore(err_act)
return
9999 continue
call psb_errpush(info,name)
call psb_error_handler(err_act)
return
contains
subroutine omega_smooth(omega,acsr)
implicit none
real(psb_spk_),intent(in) :: omega
type(psb_c_csr_sparse_mat), intent(inout) :: acsr
!
integer(psb_ipk_) :: i,j
do i=1,acsr%get_nrows()
do j=acsr%irp(i),acsr%irp(i+1)-1
if (acsr%ja(j) == i) then
acsr%val(j) = done - omega*acsr%val(j)
else
acsr%val(j) = - omega*acsr%val(j)
end if
end do
end do
end subroutine omega_smooth
end subroutine amg_c_parmatch_smth_bld
@@ -0,0 +1,160 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: amg_caggrmat_nosmth_bld.F90
!
! Subroutine: amg_caggrmat_nosmth_bld
! Version: real
!
! This routine builds a coarse-level matrix A_C from a fine-level matrix A
! by using the Galerkin approach, i.e.
!
! A_C = P_C^T A P_C,
!
! where P_C is the piecewise constant interpolation operator corresponding
! the fine-to-coarse level mapping built by amg_aggrmap_bld.
!
! The coarse-level matrix A_C is distributed among the parallel processes or
! replicated on each of them, according to the value of p%parms%coarse_mat
! specified by the user through amg_cprecinit and amg_zprecset.
! On output from this routine the entries of AC, op_prol, op_restr
! are still in "global numbering" mode; this is fixed in the calling routine
!
! For details see
! P. D'Ambra, D. di Serafino and S. Filippone, On the development of
! PSBLAS-based parallel two-level Schwarz preconditioners, Appl. Num. Math.,
! 57 (2007), 1181-1196.
!
!
! Arguments:
! a - type(psb_cspmat_type), input.
! The sparse matrix structure containing the local part of
! the fine-level matrix.
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the fine-level matrix.
! p - type(amg_c_onelev_type), input/output.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
! parms - type(amg_cml_parms), input
! Parameters controlling the choice of algorithm
! ac - type(psb_cspmat_type), output
! The coarse matrix on output
!
! ilaggr - integer, dimension(:), input
! The mapping between the row indices of the coarse-level
! matrix and the row indices of the fine-level matrix.
! ilaggr(i)=j means that node i in the adjacency graph
! of the fine-level matrix is mapped onto node j in the
! adjacency graph of the coarse-level matrix. Note that the indices
! are assumed to be shifted so as to make sure the ranges on
! the various processes do not overlap.
! nlaggr - integer, dimension(:) input
! nlaggr(i) contains the aggregates held by process i.
! op_prol - type(psb_cspmat_type), input/output
! The tentative prolongator on input, the computed prolongator on output
!
! op_restr - type(psb_cspmat_type), output
! The restrictor operator; normally, it is the transpose of the prolongator.
!
! info - integer, output.
! Error code.
!
!
subroutine amg_c_parmatch_spmm_bld(a,desc_a,ilaggr,nlaggr,parms,&
& ac,desc_ac,op_prol,op_restr,t_prol,info)
use psb_base_mod
use amg_c_inner_mod
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_spmm_bld
implicit none
! Arguments
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(inout) :: desc_a
integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:)
type(amg_cml_parms), intent(inout) :: parms
type(psb_lcspmat_type), intent(inout) :: t_prol
type(psb_cspmat_type), intent(inout) :: ac, op_prol, op_restr
type(psb_desc_type), intent(out) :: desc_ac
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_) :: err_act
type(psb_ctxt_type) :: ictxt
integer(psb_ipk_) :: np,me
character(len=20) :: name
type(psb_c_csr_sparse_mat) :: acsr
integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, nzl, ip, &
& naggr, nzt, naggrm1, naggrp1, i, k
integer(psb_ipk_) :: inaggr, nzlp
integer(psb_ipk_) :: debug_level, debug_unit
logical, parameter :: debug=.false.
name='amg_parmatch_spmm_bld'
if(psb_get_errstatus().ne.0) return
info=psb_success_
call psb_erractionsave(err_act)
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
call a%cp_to(acsr)
call amg_c_parmatch_spmm_bld_inner(acsr,desc_a,ilaggr,nlaggr,parms,&
& ac,desc_ac,op_prol,op_restr,t_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err="SPMM_BLD_INNER")
goto 9999
end if
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Done spmm_bld '
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine amg_c_parmatch_spmm_bld
@@ -0,0 +1,211 @@
!
!
! AMG4PSBLAS version 1.0
! Algebraic Multigrid Package
! based on PSBLAS (Parallel Sparse BLAS version 3.7)
!
! (C) Copyright 2021
!
! Salvatore Filippone
! Pasqua D'Ambra
! Fabio Durastante
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the AMG4PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific prior written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE AMG4PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: amg_caggrmat_nosmth_bld.F90
!
! Subroutine: amg_caggrmat_nosmth_bld
! Version: real
!
! This routine builds a coarse-level matrix A_C from a fine-level matrix A
! by using the Galerkin approach, i.e.
!
! A_C = P_C^T A P_C,
!
! where P_C is the piecewise constant interpolation operator corresponding
! the fine-to-coarse level mapping built by amg_aggrmap_bld.
!
! The coarse-level matrix A_C is distributed among the parallel processes or
! replicated on each of them, according to the value of p%parms%coarse_mat
! specified by the user through amg_cprecinit and amg_zprecset.
! On output from this routine the entries of AC, op_prol, op_restr
! are still in "global numbering" mode; this is fixed in the calling routine
!
! For details see
! P. D'Ambra, D. di Serafino and S. Filippone, On the development of
! PSBLAS-based parallel two-level Schwarz preconditioners, Appl. Num. Math.,
! 57 (2007), 1181-1196.
!
!
! Arguments:
! a - type(psb_cspmat_type), input.
! The sparse matrix structure containing the local part of
! the fine-level matrix.
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the fine-level matrix.
! p - type(amg_c_onelev_type), input/output.
! The 'one-level' data structure that will contain the local
! part of the matrix to be built as well as the information
! concerning the prolongator and its transpose.
! parms - type(amg_cml_parms), input
! Parameters controlling the choice of algorithm
! ac - type(psb_cspmat_type), output
! The coarse matrix on output
!
! ilaggr - integer, dimension(:), input
! The mapping between the row indices of the coarse-level
! matrix and the row indices of the fine-level matrix.
! ilaggr(i)=j means that node i in the adjacency graph
! of the fine-level matrix is mapped onto node j in the
! adjacency graph of the coarse-level matrix. Note that the indices
! are assumed to be shifted so as to make sure the ranges on
! the various processes do not overlap.
! nlaggr - integer, dimension(:) input
! nlaggr(i) contains the aggregates held by process i.
! op_prol - type(psb_cspmat_type), input/output
! The tentative prolongator on input, the computed prolongator on output
!
! op_restr - type(psb_cspmat_type), output
! The restrictor operator; normally, it is the transpose of the prolongator.
!
! info - integer, output.
! Error code.
!
!
subroutine amg_c_parmatch_spmm_bld_inner(a_csr,desc_a,ilaggr,nlaggr,parms,&
& ac,desc_ac,op_prol,op_restr,t_prol,info)
use psb_base_mod
use amg_c_inner_mod
use amg_c_parmatch_aggregator_mod, amg_protect_name => amg_c_parmatch_spmm_bld_inner
implicit none
! Arguments
type(psb_c_csr_sparse_mat), intent(inout) :: a_csr
type(psb_desc_type), intent(inout) :: desc_a
integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:)
type(amg_sml_parms), intent(inout) :: parms
type(psb_lcspmat_type), intent(inout) :: t_prol
type(psb_cspmat_type), intent(inout) :: ac, op_prol, op_restr
type(psb_desc_type), intent(out) :: desc_ac
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_) :: err_act
type(psb_ctxt_type) :: ictxt
integer(psb_ipk_) :: np, me, ndx
character(len=40) :: name
type(psb_lc_coo_sparse_mat) :: tmpcoo
type(psb_c_coo_sparse_mat) :: coo_prol, coo_restr
type(psb_c_csr_sparse_mat) :: ac_csr, csr_restr
type(psb_desc_type), target :: tmp_desc
type(psb_lcspmat_type) :: lac
integer(psb_ipk_) :: debug_level, debug_unit, naggr
integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, nrl, nzl, ip, &
& nzt, naggrm1, naggrp1, i, k
integer(psb_lpk_), allocatable :: ia(:),ja(:)
!integer(psb_lpk_) :: nrsave, ncsave, nzsave, nza, nrpsave, ncpsave, nzpsave
logical, parameter :: do_timings=.false., oldstyle=.false., debug=.false.
integer(psb_ipk_), save :: idx_spspmm=-1, idx_prolcnv=-1, idx_proltrans=-1, idx_asb=-1
name='amg_parmatch_spmm_bld_inner'
if(psb_get_errstatus().ne.0) return
info=psb_success_
call psb_erractionsave(err_act)
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
nglob = desc_a%get_global_rows()
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
if ((do_timings).and.(idx_spspmm==-1)) &
& idx_spspmm = psb_get_timer_idx("SPMM_BLD: spspmm ")
if ((do_timings).and.(idx_prolcnv==-1)) &
& idx_prolcnv = psb_get_timer_idx("SPMM_BLD: prolcnv ")
if ((do_timings).and.(idx_proltrans==-1)) &
& idx_proltrans = psb_get_timer_idx("SPMM_BLD: proltrans")
if ((do_timings).and.(idx_asb==-1)) &
& idx_asb = psb_get_timer_idx("SPMM_BLD: asb ")
if (do_timings) call psb_tic(idx_prolcnv)
naggr = nlaggr(me+1)
ntaggr = sum(nlaggr)
naggrm1 = sum(nlaggr(1:me))
naggrp1 = sum(nlaggr(1:me+1))
!
! Here T_PROL should be arriving with GLOBAL indices on the cols
! and LOCAL indices on the rows.
!
if (debug) write(0,*) me,' ',trim(name),' Size check on entry New: ',&
& op_prol%get_fmt(),op_prol%get_nrows(),op_prol%get_ncols(),op_prol%get_nzeros(),&
& nrow,ntaggr,naggr
call t_prol%cp_to(tmpcoo)
call psb_cdall(ictxt,desc_ac,info,nl=naggr)
nzl = tmpcoo%get_nzeros()
if (debug) write(0,*) me,' ',trim(name),' coo_prol: ',&
& tmpcoo%ia(1:min(10,nzl)),' :',tmpcoo%ja(1:min(10,nzl))
call desc_ac%indxmap%g2lip_ins(tmpcoo%ja(1:nzl),info)
call tmpcoo%set_ncols(desc_ac%get_local_cols())
call tmpcoo%cp_to_icoo(coo_prol,info)
call amg_ptap_bld(a_csr,desc_a,nlaggr,parms,ac,&
& coo_prol,desc_ac,coo_restr,info)
nzl = coo_prol%get_nzeros()
if (debug) write(0,*) me,' ',trim(name),' coo_prol: ',&
& coo_prol%ia(1:min(10,nzl)),' :',coo_prol%ja(1:min(10,nzl))
call op_prol%mv_from(coo_prol)
call op_restr%mv_from(coo_restr)
if (debug) then
write(0,*) me,' ',trim(name),' Checkpoint at exit'
call psb_barrier(ictxt)
write(0,*) me,' ',trim(name),' Checkpoint through'
end if
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,a_err='Build ac = op_restr x a3')
goto 9999
end if
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Done smooth_aggregate '
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine amg_c_parmatch_spmm_bld_inner

Some files were not shown because too many files have changed in this diff Show More