diff --git a/.VERSION b/.VERSION new file mode 100644 index 000000000..c805380c6 --- /dev/null +++ b/.VERSION @@ -0,0 +1,12 @@ +$Format:%d%n%n$ +# Fall back version, probably last release: +3.5.0 + +# PSBLAS version file. +# +# Release archive created from commit: +# $Format:%H %d$ +# $Format:Created on %ci by %cN, and$ +# $Format:signed by %GS using %GK.$ +# $Format:Signature status: %G?$ +$Format:%GG$ diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 000000000..e7ea7fd4f --- /dev/null +++ b/.gitattributes @@ -0,0 +1,47 @@ +# Handle eol correctly: Commit unix-style, checkout as native +* text=auto + +# Explicitly declare text files you want to always be normalized and converted +# to native line endings on checkout. +*.c text +*.h text +*.f90 text +*.F90 text +*.md text +*.txt text +*.sh text +*.cu text +*.x64 text + +# Denote all files that are truly binary and should not be modified. +*.mod binary +*.o binary +*.a binary +*.tar binary +*.gz binary +*.tgz binary +*.enc binary +*.pdf binary +*.png binary +*.jpg binary +*.bmp binary +*.gig binary + +# Handle windows specific files correctly +*.sln eol=crlf +*.suo eol=crlf +*.vcxproj eol=crlf +*.vcxitems eol=crlf + +# Prevent dev-ops files from making it into the release archives +.travis.yml export-ignore +.pullapprove.yml export-ignore +.gitattributes export-ignore +.gitignore export-ignore +codecov.yml export-ignore +*.enc export-ignore +.github export-ignore +.Dockerfiles export-ignore + +# Perform substitutions when `git export`ing these files +.VERSION export-subst diff --git a/.gitignore b/.gitignore index 9939c5daa..8f206ff1d 100644 --- a/.gitignore +++ b/.gitignore @@ -3,6 +3,8 @@ *.mod *~ +build + # header files generated cbind/*.h diff --git a/CMakeLists.txt b/CMakeLists.txt new file mode 100644 index 000000000..09f732f0d --- /dev/null +++ b/CMakeLists.txt @@ -0,0 +1,544 @@ +# FDEFINES : -DHAVE_LAPACK -DHAVE_FINAL -DHAVE_ISO_FORTRAN_ENV -DHAVE_FLUSH_STMT -DHAVE_VOLATILE -DSERIAL_MPI -DMPI_MOD +# CDEFINES : -DLowerUnderscore -DPtr64Bits + +#----------------------------------- +# Set oldest allowable CMake version +#----------------------------------- +cmake_minimum_required(VERSION 3.10) +cmake_policy(VERSION 3.11.1...3.13.3) + +set(CMAKE_MODULE_PATH "${CMAKE_CURRENT_LIST_DIR}/cmake") + +#---------------------------------------------- +# Define canonical CMake build types and extras +#---------------------------------------------- +set ( CMAKE_CONFIGURATION_TYPES "Debug" "Release" "MinSizeRel" "RelWithDebInfo" "CodeCoverage" ) +set ( CMAKE_BUILD_TYPE "Release" + CACHE STRING "Select which configuration to build." ) +set_property ( CACHE CMAKE_BUILD_TYPE PROPERTY STRINGS ${CMAKE_CONFIGURATION_TYPES} ) + +#----------------------------------------------------- +# Determine version from .VERSION file or git describe +#----------------------------------------------------- +include(setVersion) +set_version( + VERSION_VARIABLE PSBLAS_Version + GIT_DESCRIBE_VAR full_git_describe + CUSTOM_VERSION_FILE "${CMAKE_SOURCE_DIR}/.VERSION") +message( STATUS "Building PSBLAS version: ${full_git_describe}" ) + +#------------------------------------------ +# Name project and specify source languages +#------------------------------------------ +project(psblas + VERSION "${PSBLAS_Version}" + LANGUAGES C Fortran) + +#-------------------------------------------------- +# Set option to allow building against OpenCoarrays +#-------------------------------------------------- +if(CMAKE_Fortran_COMPILER_ID MATCHES "GNU") + option(PSBLAS_USE_OpenCoarrays "Build enabling linkage to programs using OpenCoarrays" OFF) +endif() + +#----------------------------------------------------------------- +# Define a target to create a checksummed & signed release archive +#----------------------------------------------------------------- +set(${CMAKE_PROJECT_NAME}_dist_string "${CMAKE_PROJECT_NAME}-${full_git_describe}") +if(GIT_FOUND) + add_custom_target(dist # OUTPUT "${CMAKE_BINARY_DIR}/${_package_stem_name}.tar.gz" + COMMAND "${CMAKE_COMMAND}" -P "${CMAKE_SOURCE_DIR}/cmake/makeDist.cmake" "${CMAKE_SOURCE_DIR}" "${CMAKE_BINARY_DIR}" + COMMENT "Creating source release asset, ${_package_stem_name}.tar.gz, from ${_full_git_describe} using the `git archive` command." + VERBATIM) +endif() + +#-------------------------- +# Prohibit in-source builds +#-------------------------- +include(CheckOutOfSourceBuild) + +#---------------------------------------------------- +# Define coverage flags and report untested compilers +#---------------------------------------------------- +if ("${CMAKE_Fortran_COMPILER_ID}" MATCHES "GNU" ) + set(gfortran_compiler true) + set ( CMAKE_C_FLAGS_CODECOVERAGE "-fprofile-arcs -ftest-coverage -O0" + CACHE STRING "Code coverage C compiler flags") + set ( CMAKE_Fortran_FLAGS_CODECOVERAGE "-fprofile-arcs -ftest-coverage -O0" + CACHE STRING "Code coverage Fortran compiler flags") +else() + message(WARNING + "\n" + "Attempting untested CMake build with Fortran compiler: ${CMAKE_Fortran_COMPILER_ID}. " + "Please report any failures at https://github.com/sfilippone/psblas3\n\n" + ) +endif() + +#------------------------------------ +# Fortran name mangling introspection +#------------------------------------ +include("${CMAKE_CURRENT_LIST_DIR}/cmake/CapitalizeString.cmake") +include(FortranCInterface) +CapitalizeString(${FortranCInterface_GLOBAL__CASE} fc_case) +message(STATUS "Name mangling capitalization: ${fc_case}") +message(STATUS "Name mangling fortran global suffix underscore: ${FortranCInterface_GLOBAL__SUFFIX}") +if(FortranCInterface_GLOBAL__SUFFIX STREQUAL "") + add_definitions("-D${fc_case}Case") +elseif(FortranCInterface_GLOBAL__SUFFIX STREQUAL "_") + add_definitions("-D${fc_case}Underscore") +elseif(FortranCInterface_GLOBAL__SUFFIX STREQUAL "__") + add_definitions("-D${fc_case}DoubleUnderscore") +else() + message( FATAL_ERROR "Fortran name mangling suffix, \'${FortranCInterface_GLOBAL__SUFFIX}\', unknown to PSBLAS") +endif() + +if(NOT ${WIN32}) + #---------------------------------------------- + # Determine system endian-ness and pointer size + #---------------------------------------------- + include(TestBigEndian) + TEST_BIG_ENDIAN(IS_BIG_ENDIAN) + if(IS_BIG_ENDIAN) + message( STATUS "System appears to be big endian.") + else() + message( STATUS "System appears to be little endian.") + add_definitions(-DLittleEndian) + endif() + include(CheckTypeSize) + CHECK_TYPE_SIZE("void *" VOID_P_SIZE LANGUAGE C) + if(${VOID_P_SIZE} EQUAL 8) + add_definitions(-DPtr64Bits) + endif() + message(STATUS "Have 64bit pointers") +endif() + +#----------------------------------------- +# Check for some Fortran compiler features +#----------------------------------------- +include(CheckFortranSourceCompiles) +CHECK_Fortran_SOURCE_COMPILES( + " +integer, allocatable :: a(:), b(:) +allocate(a(5)) +a = [1,2,3,4,5] +call move_alloc(from=a, to=b) +end +" + HAVE_MOVE_ALLOC + SRC_EXT f90 + ) +if(HAVE_MOVE_ALLOC) + add_definitions(-DHAVE_MOVE_ALLOC) +endif() +CHECK_Fortran_SOURCE_COMPILES( + "integer, volatile :: i ; end" + HAVE_VOLATILE + SRC_EXT f90 + ) +if(HAVE_VOLATILE) + add_definitions(-DHAVE_VOLATILE) +endif() +CHECK_Fortran_SOURCE_COMPILES( + "use ISO_FORTRAN_ENV ; end" + HAVE_ISO_FORTRAN_ENV + SRC_EXT f90 + ) +if(HAVE_ISO_FORTRAN_ENV) + add_definitions(-DHAVE_ISO_FORTRAN_ENV) +endif() +CHECK_Fortran_SOURCE_COMPILES( + "flush(5); end" + HAVE_FLUSH_STMT + SRC_EXT f90 + ) +if(HAVE_FLUSH_STMT) + add_definitions(-DHAVE_FLUSH_STMT) +endif() +CHECK_Fortran_SOURCE_COMPILES( + " +module conftest_mod + type foo + integer :: i + contains + final :: destroy_foo + end type foo + private destroy_foo +contains + subroutine destroy_foo(a) + type(foo) :: a + ! Just a test + end subroutine destroy_foo +end module conftest_mod +program conftest + use conftest_mod + type(foo) :: foovar +end program" + HAVE_FINAL + SRC_EXT f90 + ) +if(HAVE_FINAL) + add_definitions(-DHAVE_FINAL) +endif() +CHECK_Fortran_SOURCE_COMPILES( + " +program xtt + type foo + integer :: i + end type foo + type, extends(foo) :: new_foo + integer :: j + end type new_foo + class(foo), allocatable :: fooab + type(new_foo) :: nfv + integer :: info + allocate(fooab, mold=nfv, stat=info) +end program" + HAVE_MOLD + SRC_EXT f90) +if(HAVE_MOLD) + add_definitions(-DHAVE_MOLD) +endif() +CHECK_Fortran_SOURCE_COMPILES( + " +program conftest + type foo + integer :: i + end type foo + type, extends(foo) :: bar + integer j + end type bar + type(bar) :: barvar +end program " + HAVE_EXTENDS_TYPE_OF + SRC_EXT f90) +if(HAVE_EXTENDS_TYPE_OF) + add_definitions(-DHAVE_EXTENDS_TYPE_OF) +endif() +CHECK_Fortran_SOURCE_COMPILES( + " +program stt + type foo + integer :: i + end type foo + type, extends(foo) :: new_foo + integer :: j + end type new_foo + type(foo) :: foov + type(new_foo) :: nfv1, nfv2 + + + write(*,*) 'foov == nfv1? ', same_type_as(foov,nfv1) + write(*,*) 'nfv2 == nfv1? ', same_type_as(nfv2,nfv1) +end program" + HAVE_SAME_TYPE_AS + SRC_EXT f90) +if(HAVE_SAME_TYPE_AS) + add_definitions(-DHAVE_SAME_TYPE_AS) +endif() + +#---------------------------------------------------------------------------- +# Find MPI and set some flags so that FC and CC can point to gfortran and gcc +#---------------------------------------------------------------------------- +find_package( MPI ) + +if(MPI_FOUND) + #----------------------------------------------- + # Work around an issue present on fedora systems + #----------------------------------------------- + if( (MPI_C_LINK_FLAGS MATCHES "noexecstack") OR (MPI_Fortran_LINK_FLAGS MATCHES "noexecstack") ) + message ( WARNING + "The `noexecstack` linker flag was found in the MPI__LINK_FLAGS variable. This is +known to cause segmentation faults for some Fortran codes. See, e.g., +https://gcc.gnu.org/bugzilla/show_bug.cgi?id=71729 or +https://github.com/sourceryinstitute/OpenCoarrays/issues/317. + +`noexecstack` is being replaced with `execstack`" + ) + string(REPLACE "noexecstack" + "execstack" MPI_C_LINK_FLAGS_FIXED ${MPI_C_LINK_FLAGS}) + string(REPLACE "noexecstack" + "execstack" MPI_Fortran_LINK_FLAGS_FIXED ${MPI_Fortran_LINK_FLAGS}) + set(MPI_C_LINK_FLAGS "${MPI_C_LINK_FLAGS_FIXED}" CACHE STRING + "MPI C linking flags" FORCE) + set(MPI_Fortran_LINK_FLAGS "${MPI_Fortran_LINK_FLAGS_FIXED}" CACHE STRING + "MPI Fortran linking flags" FORCE) + endif() + + #---------------- + # Setup MPI flags + #---------------- + list(REMOVE_DUPLICATES MPI_Fortran_INCLUDE_PATH) + set(CMAKE_C_COMPILE_FLAGS ${CMAKE_C_COMPILE_FLAGS} ${MPI_C_COMPILE_FLAGS}) + set(CMAKE_C_LINK_FLAGS ${CMAKE_C_LINK_FLAGS} ${MPI_C_LINK_FLAGS}) + set(CMAKE_Fortran_COMPILE_FLAGS ${CMAKE_Fortran_COMPILE_FLAGS} ${MPI_Fortran_COMPILE_FLAGS}) + set(CMAKE_Fortran_LINK_FLAGS ${CMAKE_Fortran_LINK_FLAGS} ${MPI_Fortran_LINK_FLAGS}) + include_directories(BEFORE ${MPI_C_INCLUDE_PATH} ${MPI_Fortran_INCLUDE_PATH}) + + if(MPI_Fortran_HAVE_F90_MODULE OR MPI_Fortran_HAVE_F08_MODULE) + add_definitions(-DMPI_MOD) + endif() + set(SERIAL_MPI OFF) +else() + add_definitions(-DSERIAL_MPI) + add_definitions(-DMPI_MOD) + set(SERIAL_MPI ON) +endif() + +#------------------------------------------------------- +# Find and Use OpenCoarrays IFF gfortran AND options set +#------------------------------------------------------- + +if("${PSBLAS_USE_OpenCoarrays}" AND CMAKE_Fortran_COMPILER_ID MATCHES GNU) + find_package(OpenCoarrays) +endif() + +#------------------------------ +# Find Linear Algebra Libraries +#------------------------------ +if(NOT APPLE) + set(BLA_STATIC ON) +endif() +find_package(BLAS REQUIRED) +find_package(LAPACK REQUIRED) +add_definitions(-DHAVE_LAPACK) + + +#-------------------------------- +# Find METIS partitioning library +#-------------------------------- +include(${CMAKE_CURRENT_LIST_DIR}/cmake/FindMETIS.cmake) +find_package(METIS) + +#--------------------------------------------------- +# Use standardized GNU install directory conventions +#--------------------------------------------------- +include(GNUInstallDirs) +#set(mod_dir_tail "${${CMAKE_PROJECT_NAME}_dist_string}_${CMAKE_Fortran_COMPILER_ID}-${CMAKE_Fortran_COMPILER_VERSION}") +set(CMAKE_RUNTIME_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_BINDIR}/${${CMAKE_PROJECT_NAME}_dist_string}-tests") +set(CMAKE_LIBRARY_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}") +set(CMAKE_ARCHIVE_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_LIBDIR}") + +#----------------------------------- +# Turn on testing/ctest capabilities +#----------------------------------- +enable_testing() + +#------------------------------------------------------------------------------ +# Add custom properties on targets for controling number of ranks during tests +#------------------------------------------------------------------------------ +define_property(TARGET + PROPERTY MIN_RANKS + BRIEF_DOCS "Minimum allowable ranks for the test " + FULL_DOCS "Property to mark executable targets run as tests that they require at least ranks to run" + ) + +define_property(TARGET + PROPERTY POWER_2_RANKS + BRIEF_DOCS "True if test must be run with a power of 2 ranks (T/F)" + FULL_DOCS "Property to mark executable targets run as tests that they require 2^n ranks." + ) + +#----------------------------------------------------- +# Publicize installed location to other CMake projects +#----------------------------------------------------- +install(EXPORT ${CMAKE_PROJECT_NAME}-targets + DESTINATION "${CMAKE_INSTALL_LIBDIR}/cmake" +) +include(CMakePackageConfigHelpers) # standard CMake module +write_basic_package_version_file( + "${CMAKE_CURRENT_BINARY_DIR}/${CMAKE_PROJECT_NAME}ConfigVersion.cmake" + VERSION "${psblas_VERSION}" + COMPATIBILITY SameMajorVersion + ) + +configure_file("${CMAKE_SOURCE_DIR}/cmake/pkg/${CMAKE_PROJECT_NAME}Config.cmake.in" + "${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/${CMAKE_PROJECT_NAME}Config.cmake" @ONLY) + +install( + FILES + "${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/${CMAKE_PROJECT_NAME}Config.cmake" + "${CMAKE_CURRENT_BINARY_DIR}/${CMAKE_PROJECT_NAME}ConfigVersion.cmake" + DESTINATION + "${CMAKE_INSTALL_LIBDIR}/cmake/psblas" +) + +#------------------------------------------ +# Add portable unistall command to makefile +#------------------------------------------ +# Adapted from the CMake Wiki FAQ +configure_file ( "${CMAKE_SOURCE_DIR}/cmake/uninstall.cmake.in" "${CMAKE_BINARY_DIR}/uninstall.cmake" + @ONLY) + +add_custom_target ( uninstall + COMMAND ${CMAKE_COMMAND} -P "${CMAKE_BINARY_DIR}/uninstall.cmake" ) + +add_custom_target(check COMMAND ${CMAKE_CTEST_COMMAND} --output-on-failure) +# See JSON-Fortran's CMakeLists.txt file to find out how to get the check target to depend +# on the test executables + +#---------------------------------- +# Determine if we're using Open MPI +#--------------------------------- +if(MPI_FOUND) + execute_process(COMMAND ${MPIEXEC} --version + OUTPUT_VARIABLE mpi_version_out) + if (mpi_version_out MATCHES "[Oo]pen[ -][Mm][Pp][Ii]") + message( STATUS "OpenMPI detected") + set ( openmpi true ) + endif() +endif() + +#--------------------------------------- +# Add the PSBLAS libraries and utilities +#--------------------------------------- + +# Link order, left to right: +# cbind.a, util.a krylov.a prec.a base.a + +include(${CMAKE_CURRENT_LIST_DIR}/base/CMakeLists.txt) +if(WIN32) + add_library(psb_base_C STATIC ${base_source_C_files}) + target_compile_definitions(psb_base_C + PRIVATE -DWIN32 -D_LIB -DWIN64) + set_target_properties(psb_base_C + PROPERTIES + LINKER_LANGUAGE C + POSITION_INDEPENDENT_CODE TRUE) + target_link_libraries(psb_base_C + PUBLIC kernel32 user32 shell32) + add_library(base ${base_source_files}) + target_link_libraries(base + PUBLIC psb_base_C) +else() + add_library(base_C OBJECT ${base_source_C_files}) + add_library(base ${base_source_files} $) +endif() +set_target_properties(base + PROPERTIES + Fortran_MODULE_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}" + POSITION_INDEPENDENT_CODE TRUE + OUTPUT_NAME psb_base + LINKER_LANGUAGE Fortran + ) +target_include_directories(base PUBLIC + $ + $) +target_link_libraries(base + PUBLIC ${LAPACK_LINKER_FLAGS} ${LAPACK_LIBRARIES} ${LAPACK95_LIBRARIES} + PUBLIC ${BLAS_LINKER_FLAGS} ${BLAS_LIBRARIES} ${BLAS95_LIBRARIES}) + +include(${CMAKE_CURRENT_LIST_DIR}/prec/CMakeLists.txt) +add_library(prec ${prec_source_files}) +set_target_properties(prec + PROPERTIES + Fortran_MODULE_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}" + POSITION_INDEPENDENT_CODE TRUE + OUTPUT_NAME psb_prec + LINKER_LANGUAGE Fortran + ) +target_include_directories(prec PUBLIC + $ + $) +target_link_libraries(prec PUBLIC base) + +include(${CMAKE_CURRENT_LIST_DIR}/krylov/CMakeLists.txt) +add_library(krylov ${krylov_source_files}) +set_target_properties(krylov + PROPERTIES + Fortran_MODULE_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}" + POSITION_INDEPENDENT_CODE TRUE + OUTPUT_NAME psb_krylov + LINKER_LANGUAGE Fortran + ) +target_include_directories(krylov PUBLIC + $ + $) +target_link_libraries(krylov PUBLIC base prec) + +include(${CMAKE_CURRENT_LIST_DIR}/util/CMakeLists.txt) +if(WIN32) + if(METIS_FOUND) + add_library(psb_util_C STATIC ${util_source_C_files}) + target_compile_definitions(psb_util_C + PRIVATE -DWIN32 -D_LIB -DWIN64) + set_target_properties(psb_util_C + PROPERTIES + LINKER_LANGUAGE C + POSITION_INDEPENDENT_CODE TRUE) + target_link_libraries(psb_util_C + PUBLIC kernel32 user32 shell32) + endif() + add_library(util ${util_source_files}) + if(METIS_FOUND) + target_link_libraries(util + PUBLIC psb_util_C) + endif() +else() + add_library(psb_util_C OBJECT ${util_source_C_files}) + add_library(util ${util_source_files} $) +endif() +set_target_properties(util + PROPERTIES + Fortran_MODULE_DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}" + POSITION_INDEPENDENT_CODE TRUE + OUTPUT_NAME psb_util + LINKER_LANGUAGE Fortran + ) +target_include_directories(util PUBLIC + $ + $) +target_link_libraries(util PUBLIC base prec) +if(METIS_FOUND) + message(STATUS ${METIS_INCLUDES}) + target_include_directories(util + PUBLIC ${METIS_INCLUDES}) + target_include_directories(psb_util_C + PUBLIC ${METIS_INCLUDES}) + target_link_libraries(util + PUBLIC ${METIS_LIBRARIES}) + target_compile_definitions(psb_util_C + PUBLIC HAVE_METIS_) + target_compile_definitions(util + PUBLIC HAVE_METIS) +endif() + +if(MPI_FOUND) + foreach(lib base prec krylov util) + target_link_libraries(${lib} PUBLIC ${MPI_C_LIBRARIES} ${MPI_Fortran_LIBRARIES}) + endforeach() +endif() + +if(OpenCoarrays_FOUND) + foreach(lib base prec krylov util) + target_link_libraries(${lib} PUBLIC OpenCoarrays::caf_mpi_static) + endforeach() +endif() + +install(DIRECTORY "${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/" DESTINATION "${CMAKE_INSTALL_INCLUDEDIR}" + FILES_MATCHING PATTERN "*.mod") +install(TARGETS base prec krylov util + EXPORT ${CMAKE_PROJECT_NAME}-targets + DESTINATION "${CMAKE_INSTALL_LIBDIR}" + LIBRARY DESTINATION "${CMAKE_INSTALL_LIBDIR}" + ) +if(WIN32) + install(TARGETS psb_base_C + EXPORT ${CMAKE_PROJECT_NAME}-targets + DESTINATION "${CMAKE_INSTALL_LIBDIR}" + LIBRARY DESTINATION "${CMAKE_INSTALL_LIBDIR}" + ) + if(METIS_FOUND) + install(TARGETS psb_util_C + EXPORT ${CMAKE_PROJECT_NAME}-targets + DESTINATION "${CMAKE_INSTALL_LIBDIR}" + LIBRARY DESTINATION "${CMAKE_INSTALL_LIBDIR}" + ) + endif() +endif() + +#----------------- +# Add PSBLAS tests +#----------------- + +# Unit tests targeting each function, argument, and branch of code +# add_mpi_test(initialize_mpi 2 initialize_mpi) diff --git a/Changelog b/Changelog index 49cb9d172..3b9e103b0 100644 --- a/Changelog +++ b/Changelog @@ -1,5 +1,7 @@ Changelog. A lot less detailed than usual, at least for past history. +2017/12/26: Add CMake build system & rename aux directory to allow cloning on + Windows. 2017/10/02: Merged CBinding. 2017/09/30: Fixes for README, contributors, bug reporting address. 2017/08/09: New optional args to TRIL and TRIU to produce two output diff --git a/base/CMakeLists.txt b/base/CMakeLists.txt new file mode 100644 index 000000000..26e70a788 --- /dev/null +++ b/base/CMakeLists.txt @@ -0,0 +1,414 @@ +set(PSB_base_source_files +comm/internals/psi_covrl_restr.f90 +comm/internals/psi_covrl_save.f90 +comm/internals/psi_covrl_upd.f90 +comm/internals/psi_cswapdata.F90 +comm/internals/psi_cswaptran.F90 +comm/internals/psi_dovrl_restr.f90 +comm/internals/psi_dovrl_save.f90 +comm/internals/psi_dovrl_upd.f90 +comm/internals/psi_dswapdata.F90 +comm/internals/psi_dswaptran.F90 +comm/internals/psi_iovrl_restr.f90 +comm/internals/psi_iovrl_save.f90 +comm/internals/psi_iovrl_upd.f90 +comm/internals/psi_iswapdata.F90 +comm/internals/psi_iswaptran.F90 +comm/internals/psi_sovrl_restr.f90 +comm/internals/psi_sovrl_save.f90 +comm/internals/psi_sovrl_upd.f90 +comm/internals/psi_sswapdata.F90 +comm/internals/psi_sswaptran.F90 +comm/internals/psi_zovrl_restr.f90 +comm/internals/psi_zovrl_save.f90 +comm/internals/psi_zovrl_upd.f90 +comm/internals/psi_zswapdata.F90 +comm/internals/psi_zswaptran.F90 +comm/psb_cgather.f90 +comm/psb_chalo.f90 +comm/psb_covrl.f90 +comm/psb_cscatter.F90 +comm/psb_cspgather.F90 +comm/psb_dgather.f90 +comm/psb_dhalo.f90 +comm/psb_dovrl.f90 +comm/psb_dscatter.F90 +comm/psb_dspgather.F90 +comm/psb_igather.f90 +comm/psb_ihalo.f90 +comm/psb_iovrl.f90 +comm/psb_iscatter.F90 +comm/psb_sgather.f90 +comm/psb_shalo.f90 +comm/psb_sovrl.f90 +comm/psb_sscatter.F90 +comm/psb_sspgather.F90 +comm/psb_zgather.f90 +comm/psb_zhalo.f90 +comm/psb_zovrl.f90 +comm/psb_zscatter.F90 +comm/psb_zspgather.F90 +internals/psb_indx_map_fnd_owner.F90 +internals/psi_bld_tmphalo.f90 +internals/psi_bld_tmphalo.f90 +internals/psi_bld_tmpovrl.f90 +internals/psi_bld_tmpovrl.f90 +internals/psi_compute_size.f90 +internals/psi_compute_size.f90 +internals/psi_crea_bnd_elem.f90 +internals/psi_crea_bnd_elem.f90 +internals/psi_crea_index.f90 +internals/psi_crea_index.f90 +internals/psi_crea_ovr_elem.f90 +internals/psi_crea_ovr_elem.f90 +internals/psi_desc_impl.f90 +internals/psi_desc_impl.f90 +internals/psi_desc_index.F90 +internals/psi_dl_check.f90 +internals/psi_dl_check.f90 +internals/psi_exist_ovr_elem.f90 +internals/psi_exist_ovr_elem.f90 +internals/psi_extrct_dl.F90 +internals/psi_fnd_owner.F90 +internals/psi_list_search.f90 +internals/psi_list_search.f90 +internals/psi_sort_dl.f90 +internals/psi_sort_dl.f90 +internals/psi_srtlist.f90 +internals/psi_srtlist.f90 +modules/auxil/psb_c_sort_mod.f90 +modules/auxil/psb_d_sort_mod.f90 +modules/auxil/psb_hash_mod.f90 +modules/auxil/psb_i_sort_mod.f90 +modules/auxil/psb_ip_reord_mod.f90 +modules/auxil/psb_s_sort_mod.f90 +modules/auxil/psb_sort_mod.f90 +modules/auxil/psb_string_mod.f90 +modules/auxil/psb_z_sort_mod.f90 +modules/auxil/psi_c_serial_mod.f90 +modules/auxil/psi_d_serial_mod.f90 +modules/auxil/psi_i_serial_mod.f90 +modules/auxil/psi_s_serial_mod.f90 +modules/auxil/psi_serial_mod.f90 +modules/auxil/psi_z_serial_mod.f90 +modules/comm/psb_base_linmap_mod.f90 +modules/comm/psb_c_comm_mod.f90 +modules/comm/psb_c_linmap_mod.f90 +modules/comm/psb_comm_mod.f90 +modules/comm/psb_d_comm_mod.f90 +modules/comm/psb_d_linmap_mod.f90 +modules/comm/psb_i_comm_mod.f90 +modules/comm/psb_linmap_mod.f90 +modules/comm/psb_s_comm_mod.f90 +modules/comm/psb_s_linmap_mod.f90 +modules/comm/psb_z_comm_mod.f90 +modules/comm/psb_z_linmap_mod.f90 +modules/desc/psb_desc_const_mod.f90 +modules/desc/psb_desc_mod.F90 +modules/desc/psb_gen_block_map_mod.f90 +modules/desc/psb_glist_map_mod.f90 +modules/desc/psb_hash_map_mod.f90 +modules/desc/psb_indx_map_mod.f90 +modules/desc/psb_list_map_mod.f90 +modules/desc/psb_repl_map_mod.f90 +modules/error.f90 +modules/psb_base_mod.f90 +modules/psb_check_mod.f90 +modules/psb_const_mod.F90 +modules/psb_error_impl.F90 +modules/psb_error_mod.F90 +modules/psb_penv_mod.F90 +modules/psb_realloc_mod.F90 +modules/psblas/psb_c_psblas_mod.F90 +modules/psblas/psb_d_psblas_mod.F90 +modules/psblas/psb_psblas_mod.f90 +modules/psblas/psb_s_psblas_mod.F90 +modules/psblas/psb_z_psblas_mod.F90 +modules/psi_bcast_mod.F90 +modules/psi_c_mod.f90 +modules/psi_comm_buffers_mod.F90 +modules/psi_d_mod.f90 +modules/psi_i_mod.f90 +modules/psi_mod.f90 +modules/psi_p2p_mod.F90 +modules/psi_penv_mod.F90 +modules/psi_reduce_mod.F90 +modules/psi_s_mod.f90 +modules/psi_z_mod.f90 +modules/serial/psb_base_mat_mod.f90 +modules/serial/psb_c_base_mat_mod.f90 +modules/serial/psb_c_base_vect_mod.f90 +modules/serial/psb_c_csc_mat_mod.f90 +modules/serial/psb_c_csr_mat_mod.f90 +modules/serial/psb_c_mat_mod.f90 +modules/serial/psb_c_serial_mod.f90 +modules/serial/psb_c_vect_mod.F90 +modules/serial/psb_d_base_mat_mod.f90 +modules/serial/psb_d_base_vect_mod.f90 +modules/serial/psb_d_csc_mat_mod.f90 +modules/serial/psb_d_csr_mat_mod.f90 +modules/serial/psb_d_mat_mod.f90 +modules/serial/psb_d_serial_mod.f90 +modules/serial/psb_d_vect_mod.F90 +modules/serial/psb_i_base_vect_mod.f90 +modules/serial/psb_i_vect_mod.F90 +modules/serial/psb_mat_mod.f90 +modules/serial/psb_s_base_mat_mod.f90 +modules/serial/psb_s_base_vect_mod.f90 +modules/serial/psb_s_csc_mat_mod.f90 +modules/serial/psb_s_csr_mat_mod.f90 +modules/serial/psb_s_mat_mod.f90 +modules/serial/psb_s_serial_mod.f90 +modules/serial/psb_s_vect_mod.F90 +modules/serial/psb_serial_mod.f90 +modules/serial/psb_vect_mod.f90 +modules/serial/psb_z_base_mat_mod.f90 +modules/serial/psb_z_base_vect_mod.f90 +modules/serial/psb_z_csc_mat_mod.f90 +modules/serial/psb_z_csr_mat_mod.f90 +modules/serial/psb_z_mat_mod.f90 +modules/serial/psb_z_serial_mod.f90 +modules/serial/psb_z_vect_mod.F90 +modules/tools/psb_c_tools_mod.f90 +modules/tools/psb_cd_tools_mod.f90 +modules/tools/psb_d_tools_mod.f90 +modules/tools/psb_i_tools_mod.f90 +modules/tools/psb_s_tools_mod.f90 +modules/tools/psb_tools_mod.f90 +modules/tools/psb_z_tools_mod.f90 +psblas/psb_camax.f90 +psblas/psb_casum.f90 +psblas/psb_caxpby.f90 +psblas/psb_cdot.f90 +psblas/psb_cnrm2.f90 +psblas/psb_cnrmi.f90 +psblas/psb_cspmm.f90 +psblas/psb_cspnrm1.f90 +psblas/psb_cspsm.f90 +psblas/psb_damax.f90 +psblas/psb_dasum.f90 +psblas/psb_daxpby.f90 +psblas/psb_ddot.f90 +psblas/psb_dnrm2.f90 +psblas/psb_dnrmi.f90 +psblas/psb_dspmm.f90 +psblas/psb_dspnrm1.f90 +psblas/psb_dspsm.f90 +psblas/psb_samax.f90 +psblas/psb_sasum.f90 +psblas/psb_saxpby.f90 +psblas/psb_sdot.f90 +psblas/psb_snrm2.f90 +psblas/psb_snrmi.f90 +psblas/psb_sspmm.f90 +psblas/psb_sspnrm1.f90 +psblas/psb_sspsm.f90 +psblas/psb_zamax.f90 +psblas/psb_zasum.f90 +psblas/psb_zaxpby.f90 +psblas/psb_zdot.f90 +psblas/psb_znrm2.f90 +psblas/psb_znrmi.f90 +psblas/psb_zspmm.f90 +psblas/psb_zspnrm1.f90 +psblas/psb_zspsm.f90 +serial/impl/psb_s_coo_impl.f90 +serial/impl/psb_s_csr_impl.f90 +serial/impl/psb_c_base_mat_impl.F90 +serial/impl/psb_base_mat_impl.f90 +serial/impl/psb_d_mat_impl.F90 +serial/impl/psb_d_csc_impl.f90 +serial/impl/psb_c_mat_impl.F90 +serial/impl/psb_c_csc_impl.f90 +serial/impl/psb_z_csc_impl.f90 +serial/impl/psb_z_mat_impl.F90 +serial/impl/psb_d_base_mat_impl.F90 +serial/impl/psb_z_base_mat_impl.F90 +serial/impl/psb_d_coo_impl.f90 +serial/impl/psb_c_coo_impl.f90 +serial/impl/psb_s_base_mat_impl.F90 +serial/impl/psb_z_csr_impl.f90 +serial/impl/psb_d_csr_impl.f90 +serial/impl/psb_c_csr_impl.f90 +serial/impl/psb_z_coo_impl.f90 +serial/impl/psb_s_mat_impl.F90 +serial/impl/psb_s_csc_impl.f90 +serial/psi_z_serial_impl.f90 +serial/psb_sasum_s.f90 +serial/psb_damax_s.f90 +serial/psb_dnumbmm.f90 +serial/psb_drwextd.f90 +serial/psb_zgeprt.f90 +serial/psb_dgelp.f90 +serial/psb_csymbmm.f90 +serial/psb_zamax_s.f90 +serial/psb_sspspmm.f90 +serial/psb_dgeprt.f90 +serial/psb_spdot_srtd.f90 +serial/psb_znumbmm.f90 +serial/psb_zrwextd.f90 +serial/psb_lsame.f90 +serial/sort/psb_s_qsort_impl.f90 +serial/sort/psi_lcx_mod.f90 +serial/sort/psb_z_hsort_impl.f90 +serial/sort/psb_i_msort_impl.f90 +serial/sort/psb_s_isort_impl.f90 +serial/sort/psi_acx_mod.f90 +serial/sort/psb_z_isort_impl.f90 +serial/sort/psb_z_qsort_impl.f90 +serial/sort/psb_s_hsort_impl.f90 +serial/sort/psb_d_msort_impl.f90 +serial/sort/psb_c_msort_impl.f90 +serial/sort/psb_d_hsort_impl.f90 +serial/sort/psb_c_hsort_impl.f90 +serial/sort/psb_d_isort_impl.f90 +serial/sort/psb_c_qsort_impl.f90 +serial/sort/psb_c_isort_impl.f90 +serial/sort/psb_d_qsort_impl.f90 +serial/sort/psb_i_isort_impl.f90 +serial/sort/psi_alcx_mod.f90 +serial/sort/psb_i_qsort_impl.f90 +serial/sort/psb_s_msort_impl.f90 +serial/sort/psb_i_hsort_impl.f90 +serial/sort/psb_z_msort_impl.f90 +serial/psb_cgeprt.f90 +serial/smmp.f90 +serial/psb_zsymbmm.f90 +serial/psi_s_serial_impl.f90 +serial/psi_i_serial_impl.f90 +serial/psb_zgelp.f90 +serial/psb_dsymbmm.f90 +serial/psb_aspxpby.f90 +serial/psb_crwextd.f90 +serial/psb_cnumbmm.f90 +serial/psb_sgelp.f90 +serial/psb_camax_s.f90 +serial/psb_cgelp.f90 +serial/psi_d_serial_impl.f90 +serial/psb_dasum_s.f90 +serial/psb_zspspmm.f90 +serial/psb_samax_s.f90 +serial/psb_snumbmm.f90 +serial/psb_srwextd.f90 +serial/psb_dspspmm.f90 +serial/psi_c_serial_impl.f90 +serial/psb_sgeprt.f90 +serial/psb_zasum_s.f90 +serial/psb_cspspmm.f90 +serial/psb_spge_dot.f90 +serial/psb_casum_s.f90 +serial/psb_ssymbmm.f90 +tools/psb_c_map.f90 +tools/psb_c_map.f90 +tools/psb_callc.f90 +tools/psb_casb.f90 +tools/psb_ccdbldext.F90 +tools/psb_cd_inloc.f90 +tools/psb_cd_lstext.f90 +tools/psb_cd_lstext.f90 +tools/psb_cd_reinit.f90 +tools/psb_cd_set_bld.f90 +tools/psb_cd_switch_ovl_indxmap.f90 +tools/psb_cd_switch_ovl_indxmap.f90 +tools/psb_cdall.f90 +tools/psb_cdals.f90 +tools/psb_cdals.f90 +tools/psb_cdalv.f90 +tools/psb_cdalv.f90 +tools/psb_cdcpy.F90 +tools/psb_cdins.f90 +tools/psb_cdins.f90 +tools/psb_cdprt.f90 +tools/psb_cdren.f90 +tools/psb_cdrep.f90 +tools/psb_cdrep.f90 +tools/psb_cfree.f90 +tools/psb_cins.f90 +tools/psb_cins.f90 +tools/psb_cspalloc.f90 +tools/psb_cspasb.f90 +tools/psb_cspasb.f90 +tools/psb_cspfree.f90 +tools/psb_cspfree.f90 +tools/psb_csphalo.F90 +tools/psb_cspins.f90 +tools/psb_cspins.f90 +tools/psb_csprn.f90 +tools/psb_csprn.f90 +tools/psb_d_map.f90 +tools/psb_dallc.f90 +tools/psb_dallc.f90 +tools/psb_dasb.f90 +tools/psb_dasb.f90 +tools/psb_dcdbldext.F90 +tools/psb_dfree.f90 +tools/psb_dfree.f90 +tools/psb_dins.f90 +tools/psb_dins.f90 +tools/psb_dspalloc.f90 +tools/psb_dspalloc.f90 +tools/psb_dspasb.f90 +tools/psb_dspfree.f90 +tools/psb_dspfree.f90 +tools/psb_dsphalo.F90 +tools/psb_dspins.f90 +tools/psb_dsprn.f90 +tools/psb_dsprn.f90 +tools/psb_get_overlap.f90 +tools/psb_get_overlap.f90 +tools/psb_glob_to_loc.f90 +tools/psb_iallc.f90 +tools/psb_iallc.f90 +tools/psb_iasb.f90 +tools/psb_iasb.f90 +tools/psb_icdasb.F90 +tools/psb_ifree.f90 +tools/psb_ifree.f90 +tools/psb_iins.f90 +tools/psb_loc_to_glob.f90 +tools/psb_s_map.f90 +tools/psb_sallc.f90 +tools/psb_sasb.f90 +tools/psb_sasb.f90 +tools/psb_scdbldext.F90 +tools/psb_sfree.f90 +tools/psb_sins.f90 +tools/psb_sins.f90 +tools/psb_sspalloc.f90 +tools/psb_sspasb.f90 +tools/psb_sspasb.f90 +tools/psb_sspfree.f90 +tools/psb_sspfree.f90 +tools/psb_ssphalo.F90 +tools/psb_sspins.f90 +tools/psb_ssprn.f90 +tools/psb_z_map.f90 +tools/psb_z_map.f90 +tools/psb_zallc.f90 +tools/psb_zallc.f90 +tools/psb_zasb.f90 +tools/psb_zcdbldext.F90 +tools/psb_zfree.f90 +tools/psb_zfree.f90 +tools/psb_zins.f90 +tools/psb_zins.f90 +tools/psb_zspalloc.f90 +tools/psb_zspasb.f90 +tools/psb_zspfree.f90 +tools/psb_zspfree.f90 +tools/psb_zsphalo.F90 +tools/psb_zspins.f90 +tools/psb_zspins.f90 +tools/psb_zsprn.f90 +) +foreach(file IN LISTS PSB_base_source_files) + list(APPEND base_source_files ${CMAKE_CURRENT_LIST_DIR}/${file}) +endforeach() + +list(APPEND PSB_base_source_C_files modules/cutil.c) +if (SERIAL_MPI) + list(APPEND PSB_base_source_C_files modules/fakempi.c) +endif() +foreach(file IN LISTS PSB_base_source_C_files) + list(APPEND base_source_C_files ${CMAKE_CURRENT_LIST_DIR}/${file}) +endforeach() diff --git a/base/comm/internals/psi_dswapdata.F90 b/base/comm/internals/psi_dswapdata.F90 index f00d21b94..9b56e2824 100644 --- a/base/comm/internals/psi_dswapdata.F90 +++ b/base/comm/internals/psi_dswapdata.F90 @@ -365,8 +365,9 @@ subroutine psi_dswapidxm(iictxt,iicomm,flag,n,beta,y,idx, & ! Then I post all the blocking sends +#ifndef SERIAL_MPI if (usersend) call mpi_barrier(icomm,iret) - +#endif pnti = 1 snd_pt = 1 rcv_pt = 1 @@ -376,17 +377,20 @@ subroutine psi_dswapidxm(iictxt,iicomm,flag,n,beta,y,idx, & nesd = idx(pnti+nerv+psb_n_elem_send_) p2ptag = psb_double_swap_tag - if ((nesd>0).and.(proc_to_comm /= me)) then - if (usersend) then + if ((nesd>0).and.(proc_to_comm /= me)) then +#ifndef SERIAL_MPI + if (usersend) then call mpi_rsend(sndbuf(snd_pt),n*nesd,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag,icomm,iret) else +#endif call mpi_send(sndbuf(snd_pt),n*nesd,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag,icomm,iret) +#ifndef SERIAL_MPI end if - +#endif if(iret /= mpi_success) then ierr(1) = iret info=psb_err_mpi_error_ @@ -855,8 +859,9 @@ subroutine psi_dswapidxv(iictxt,iicomm,flag,beta,y,idx, & ! Then I post all the blocking sends +#ifndef SERIAL_MPI if (usersend) call mpi_barrier(icomm,iret) - +#endif pnti = 1 snd_pt = 1 rcv_pt = 1 @@ -868,16 +873,19 @@ subroutine psi_dswapidxv(iictxt,iicomm,flag,beta,y,idx, & p2ptag = psb_double_swap_tag if ((nesd>0).and.(proc_to_comm /= me)) then +#ifndef SERIAL_MPI if (usersend) then call mpi_rsend(sndbuf(snd_pt),nesd,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag,icomm,iret) else +#endif call mpi_send(sndbuf(snd_pt),nesd,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag,icomm,iret) +#ifndef SERIAL_MPI end if - +#endif if(iret /= mpi_success) then ierr(1) = iret info=psb_err_mpi_error_ diff --git a/base/comm/internals/psi_dswaptran.F90 b/base/comm/internals/psi_dswaptran.F90 index 604b345e6..2b328a2eb 100644 --- a/base/comm/internals/psi_dswaptran.F90 +++ b/base/comm/internals/psi_dswaptran.F90 @@ -376,8 +376,9 @@ subroutine psi_dtranidxm(iictxt,iicomm,flag,n,beta,y,idx,& ! Then I post all the blocking sends +#ifndef SERIAL_MPI if (usersend) call mpi_barrier(icomm,iret) - +#endif pnti = 1 snd_pt = 1 rcv_pt = 1 @@ -388,16 +389,19 @@ subroutine psi_dtranidxm(iictxt,iicomm,flag,n,beta,y,idx,& if ((nerv>0).and.(proc_to_comm /= me)) then p2ptag = psb_double_swap_tag +#ifndef SERIAL_MPI if (usersend) then call mpi_rsend(rcvbuf(rcv_pt),n*nerv,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag,icomm,iret) else +#endif call mpi_send(rcvbuf(rcv_pt),n*nerv,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag,icomm,iret) +#ifndef SERIAL_MPI end if - +#endif if(iret /= mpi_success) then ierr(1) = iret info=psb_err_mpi_error_ @@ -874,8 +878,9 @@ subroutine psi_dtranidxv(iictxt,iicomm,flag,beta,y,idx,& ! Then I post all the blocking sends +#ifndef SERIAL_MPI if (usersend) call mpi_barrier(icomm,iret) - +#endif pnti = 1 snd_pt = 1 rcv_pt = 1 @@ -886,16 +891,19 @@ subroutine psi_dtranidxv(iictxt,iicomm,flag,beta,y,idx,& if ((nerv>0).and.(proc_to_comm /= me)) then p2ptag = psb_double_swap_tag +#ifndef SERIAL_MPI if (usersend) then call mpi_rsend(rcvbuf(rcv_pt),nerv,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag, icomm,iret) else +#endif call mpi_send(rcvbuf(rcv_pt),nerv,& & psb_mpi_r_dpk_,prcid(i),& & p2ptag, icomm,iret) +#ifndef SERIAL_MPI end if - +#endif if(iret /= mpi_success) then ierr(1) = iret info=psb_err_mpi_error_ diff --git a/base/modules/Makefile b/base/modules/Makefile index ec53f942c..cde7ee900 100644 --- a/base/modules/Makefile +++ b/base/modules/Makefile @@ -2,10 +2,10 @@ include ../../Make.inc BASIC_MODS= psb_const_mod.o psb_error_mod.o psb_realloc_mod.o COMMINT=psi_comm_buffers_mod.o psi_penv_mod.o psi_bcast_mod.o psi_reduce_mod.o psi_p2p_mod.o -UTIL_MODS = aux/psb_string_mod.o desc/psb_desc_const_mod.o desc/psb_indx_map_mod.o\ +UTIL_MODS = auxil/psb_string_mod.o desc/psb_desc_const_mod.o desc/psb_indx_map_mod.o\ desc/psb_gen_block_map_mod.o desc/psb_list_map_mod.o desc/psb_repl_map_mod.o\ desc/psb_glist_map_mod.o desc/psb_hash_map_mod.o \ - desc/psb_desc_mod.o aux/psb_sort_mod.o \ + desc/psb_desc_mod.o auxil/psb_sort_mod.o \ serial/psb_s_serial_mod.o serial/psb_d_serial_mod.o serial/psb_c_serial_mod.o serial/psb_z_serial_mod.o \ serial/psb_serial_mod.o \ tools/psb_cd_tools_mod.o tools/psb_i_tools_mod.o tools/psb_s_tools_mod.o tools/psb_d_tools_mod.o\ @@ -23,13 +23,13 @@ UTIL_MODS = aux/psb_string_mod.o desc/psb_desc_const_mod.o desc/psb_indx_map_mod serial/psb_vect_mod.o\ psblas/psb_s_psblas_mod.o psblas/psb_c_psblas_mod.o \ psblas/psb_d_psblas_mod.o psblas/psb_z_psblas_mod.o psblas/psb_psblas_mod.o \ - aux/psi_serial_mod.o aux/psi_i_serial_mod.o \ - aux/psi_s_serial_mod.o aux/psi_d_serial_mod.o aux/psi_c_serial_mod.o aux/psi_z_serial_mod.o \ + auxil/psi_serial_mod.o auxil/psi_i_serial_mod.o \ + auxil/psi_s_serial_mod.o auxil/psi_d_serial_mod.o auxil/psi_c_serial_mod.o auxil/psi_z_serial_mod.o \ psi_mod.o psi_i_mod.o psi_s_mod.o psi_d_mod.o psi_c_mod.o psi_z_mod.o\ - aux/psb_ip_reord_mod.o\ - aux/psb_i_sort_mod.o aux/psb_s_sort_mod.o aux/psb_d_sort_mod.o \ - aux/psb_c_sort_mod.o aux/psb_z_sort_mod.o \ - psb_check_mod.o aux/psb_hash_mod.o\ + auxil/psb_ip_reord_mod.o\ + auxil/psb_i_sort_mod.o auxil/psb_s_sort_mod.o auxil/psb_d_sort_mod.o \ + auxil/psb_c_sort_mod.o auxil/psb_z_sort_mod.o \ + psb_check_mod.o auxil/psb_hash_mod.o\ serial/psb_base_mat_mod.o serial/psb_mat_mod.o\ serial/psb_s_base_mat_mod.o serial/psb_s_csr_mat_mod.o serial/psb_s_csc_mat_mod.o serial/psb_s_mat_mod.o \ serial/psb_d_base_mat_mod.o serial/psb_d_csr_mat_mod.o serial/psb_d_csc_mat_mod.o serial/psb_d_mat_mod.o \ @@ -60,24 +60,24 @@ psi_penv_mod.o: psi_comm_buffers_mod.o psi_bcast_mod.o psi_reduce_mod.o psi_p2p_mod.o: psi_penv_mod.o -aux/psb_string_mod.o desc/psb_desc_const_mod.o psi_comm_buffers_mod.o: psb_const_mod.o -aux/psb_hash_mod.o: psb_realloc_mod.o psb_const_mod.o -aux/psb_i_sort_mod.o aux/psb_s_sort_mod.o aux/psb_d_sort_mod.o aux/psb_c_sort_mod.o aux/psb_z_sort_mod.o \ -aux/psb_ip_reord_mod.o aux/psi_serial_mod.o aux/psb_sort_mod.o: $(BASIC_MODS) -aux/psb_sort_mod.o: aux/psb_i_sort_mod.o aux/psb_s_sort_mod.o aux/psb_d_sort_mod.o \ - aux/psb_c_sort_mod.o aux/psb_z_sort_mod.o aux/psb_ip_reord_mod.o aux/psi_serial_mod.o -aux/psi_serial_mod.o: aux/psi_i_serial_mod.o \ - aux/psi_s_serial_mod.o aux/psi_d_serial_mod.o aux/psi_c_serial_mod.o aux/psi_z_serial_mod.o -aux/psi_i_serial_mod.o aux/psi_s_serial_mod.o aux/psi_d_serial_mod.o aux/psi_c_serial_mod.o aux/psi_z_serial_mod.o: psb_const_mod.o +auxil/psb_string_mod.o desc/psb_desc_const_mod.o psi_comm_buffers_mod.o: psb_const_mod.o +auxil/psb_hash_mod.o: psb_realloc_mod.o psb_const_mod.o +auxil/psb_i_sort_mod.o auxil/psb_s_sort_mod.o auxil/psb_d_sort_mod.o auxil/psb_c_sort_mod.o auxil/psb_z_sort_mod.o \ +auxil/psb_ip_reord_mod.o auxil/psi_serial_mod.o auxil/psb_sort_mod.o: $(BASIC_MODS) +auxil/psb_sort_mod.o: auxil/psb_i_sort_mod.o auxil/psb_s_sort_mod.o auxil/psb_d_sort_mod.o \ + auxil/psb_c_sort_mod.o auxil/psb_z_sort_mod.o auxil/psb_ip_reord_mod.o auxil/psi_serial_mod.o +auxil/psi_serial_mod.o: auxil/psi_i_serial_mod.o \ + auxil/psi_s_serial_mod.o auxil/psi_d_serial_mod.o auxil/psi_c_serial_mod.o auxil/psi_z_serial_mod.o +auxil/psi_i_serial_mod.o auxil/psi_s_serial_mod.o auxil/psi_d_serial_mod.o auxil/psi_c_serial_mod.o auxil/psi_z_serial_mod.o: psb_const_mod.o -serial/psb_base_mat_mod.o: aux/psi_serial_mod.o +serial/psb_base_mat_mod.o: auxil/psi_serial_mod.o serial/psb_s_base_mat_mod.o serial/psb_d_base_mat_mod.o serial/psb_c_base_mat_mod.o serial/psb_z_base_mat_mod.o: serial/psb_base_mat_mod.o serial/psb_s_base_mat_mod.o: serial/psb_s_base_vect_mod.o serial/psb_d_base_mat_mod.o: serial/psb_d_base_vect_mod.o serial/psb_c_base_mat_mod.o: serial/psb_c_base_vect_mod.o serial/psb_z_base_mat_mod.o: serial/psb_z_base_vect_mod.o serial/psb_c_base_vect_mod.o serial/psb_s_base_vect_mod.o serial/psb_d_base_vect_mod.o serial/psb_z_base_vect_mod.o: serial/psb_i_base_vect_mod.o -serial/psb_i_base_vect_mod.o serial/psb_c_base_vect_mod.o serial/psb_s_base_vect_mod.o serial/psb_d_base_vect_mod.o serial/psb_z_base_vect_mod.o: aux/psi_serial_mod.o psb_realloc_mod.o +serial/psb_i_base_vect_mod.o serial/psb_c_base_vect_mod.o serial/psb_s_base_vect_mod.o serial/psb_d_base_vect_mod.o serial/psb_z_base_vect_mod.o: auxil/psi_serial_mod.o psb_realloc_mod.o serial/psb_s_mat_mod.o: serial/psb_s_base_mat_mod.o serial/psb_s_csr_mat_mod.o serial/psb_s_csc_mat_mod.o serial/psb_s_vect_mod.o serial/psb_d_mat_mod.o: serial/psb_d_base_mat_mod.o serial/psb_d_csr_mat_mod.o serial/psb_d_csc_mat_mod.o serial/psb_d_vect_mod.o serial/psb_i_vect_mod.o serial/psb_c_mat_mod.o: serial/psb_c_base_mat_mod.o serial/psb_c_csr_mat_mod.o serial/psb_c_csc_mat_mod.o serial/psb_c_vect_mod.o @@ -93,15 +93,15 @@ serial/psb_s_vect_mod.o: serial/psb_s_base_vect_mod.o serial/psb_i_vect_mod.o serial/psb_d_vect_mod.o: serial/psb_d_base_vect_mod.o serial/psb_i_vect_mod.o serial/psb_c_vect_mod.o: serial/psb_c_base_vect_mod.o serial/psb_i_vect_mod.o serial/psb_z_vect_mod.o: serial/psb_z_base_vect_mod.o serial/psb_i_vect_mod.o -serial/psb_s_serial_mod.o serial/psb_d_serial_mod.o serial/psb_c_serial_mod.o serial/psb_z_serial_mod.o: serial/psb_mat_mod.o aux/psb_string_mod.o aux/psb_sort_mod.o aux/psi_serial_mod.o +serial/psb_s_serial_mod.o serial/psb_d_serial_mod.o serial/psb_c_serial_mod.o serial/psb_z_serial_mod.o: serial/psb_mat_mod.o auxil/psb_string_mod.o auxil/psb_sort_mod.o auxil/psi_serial_mod.o serial/psb_vect_mod.o: serial/psb_i_vect_mod.o serial/psb_d_vect_mod.o serial/psb_s_vect_mod.o serial/psb_c_vect_mod.o serial/psb_z_vect_mod.o error.o psb_realloc_mod.o: psb_error_mod.o psb_error_impl.o: psb_penv_mod.o -psb_spmat_type.o: aux/psb_string_mod.o aux/psb_sort_mod.o +psb_spmat_type.o: auxil/psb_string_mod.o auxil/psb_sort_mod.o desc/psb_desc_mod.o: psb_penv_mod.o psb_realloc_mod.o\ - aux/psb_hash_mod.o desc/psb_hash_map_mod.o desc/psb_list_map_mod.o \ + auxil/psb_hash_mod.o desc/psb_hash_map_mod.o desc/psb_list_map_mod.o \ desc/psb_repl_map_mod.o desc/psb_gen_block_map_mod.o desc/psb_desc_const_mod.o\ desc/psb_indx_map_mod.o serial/psb_i_vect_mod.o psi_i_mod.o: desc/psb_desc_mod.o serial/psb_i_vect_mod.o @@ -109,16 +109,16 @@ psi_s_mod.o: desc/psb_desc_mod.o serial/psb_s_vect_mod.o psi_d_mod.o: desc/psb_desc_mod.o serial/psb_d_vect_mod.o psi_c_mod.o: desc/psb_desc_mod.o serial/psb_c_vect_mod.o psi_z_mod.o: desc/psb_desc_mod.o serial/psb_z_vect_mod.o -psi_mod.o: psb_penv_mod.o desc/psb_desc_mod.o aux/psi_serial_mod.o serial/psb_serial_mod.o\ +psi_mod.o: psb_penv_mod.o desc/psb_desc_mod.o auxil/psi_serial_mod.o serial/psb_serial_mod.o\ psi_i_mod.o psi_s_mod.o psi_d_mod.o psi_c_mod.o psi_z_mod.o desc/psb_indx_map_mod.o: desc/psb_desc_const_mod.o psb_error_mod.o psb_penv_mod.o desc/psb_hash_map_mod.o desc/psb_list_map_mod.o desc/psb_repl_map_mod.o desc/psb_gen_block_map_mod.o:\ desc/psb_indx_map_mod.o desc/psb_desc_const_mod.o \ - aux/psb_sort_mod.o psb_penv_mod.o + auxil/psb_sort_mod.o psb_penv_mod.o desc/psb_glist_map_mod.o: desc/psb_list_map_mod.o -desc/psb_hash_map_mod.o: aux/psb_hash_mod.o aux/psb_sort_mod.o -desc/psb_gen_block_map_mod.o: aux/psb_hash_mod.o +desc/psb_hash_map_mod.o: auxil/psb_hash_mod.o auxil/psb_sort_mod.o +desc/psb_gen_block_map_mod.o: auxil/psb_hash_mod.o psb_check_mod.o: desc/psb_desc_mod.o diff --git a/base/modules/aux/psb_c_sort_mod.f90 b/base/modules/auxil/psb_c_sort_mod.f90 similarity index 100% rename from base/modules/aux/psb_c_sort_mod.f90 rename to base/modules/auxil/psb_c_sort_mod.f90 diff --git a/base/modules/aux/psb_d_sort_mod.f90 b/base/modules/auxil/psb_d_sort_mod.f90 similarity index 100% rename from base/modules/aux/psb_d_sort_mod.f90 rename to base/modules/auxil/psb_d_sort_mod.f90 diff --git a/base/modules/aux/psb_hash_mod.f90 b/base/modules/auxil/psb_hash_mod.f90 similarity index 100% rename from base/modules/aux/psb_hash_mod.f90 rename to base/modules/auxil/psb_hash_mod.f90 diff --git a/base/modules/aux/psb_i_sort_mod.f90 b/base/modules/auxil/psb_i_sort_mod.f90 similarity index 100% rename from base/modules/aux/psb_i_sort_mod.f90 rename to base/modules/auxil/psb_i_sort_mod.f90 diff --git a/base/modules/aux/psb_ip_reord_mod.f90 b/base/modules/auxil/psb_ip_reord_mod.f90 similarity index 100% rename from base/modules/aux/psb_ip_reord_mod.f90 rename to base/modules/auxil/psb_ip_reord_mod.f90 diff --git a/base/modules/aux/psb_s_sort_mod.f90 b/base/modules/auxil/psb_s_sort_mod.f90 similarity index 100% rename from base/modules/aux/psb_s_sort_mod.f90 rename to base/modules/auxil/psb_s_sort_mod.f90 diff --git a/base/modules/aux/psb_sort_mod.f90 b/base/modules/auxil/psb_sort_mod.f90 similarity index 100% rename from base/modules/aux/psb_sort_mod.f90 rename to base/modules/auxil/psb_sort_mod.f90 diff --git a/base/modules/aux/psb_string_mod.f90 b/base/modules/auxil/psb_string_mod.f90 similarity index 100% rename from base/modules/aux/psb_string_mod.f90 rename to base/modules/auxil/psb_string_mod.f90 diff --git a/base/modules/aux/psb_z_sort_mod.f90 b/base/modules/auxil/psb_z_sort_mod.f90 similarity index 100% rename from base/modules/aux/psb_z_sort_mod.f90 rename to base/modules/auxil/psb_z_sort_mod.f90 diff --git a/base/modules/aux/psi_c_serial_mod.f90 b/base/modules/auxil/psi_c_serial_mod.f90 similarity index 100% rename from base/modules/aux/psi_c_serial_mod.f90 rename to base/modules/auxil/psi_c_serial_mod.f90 diff --git a/base/modules/aux/psi_d_serial_mod.f90 b/base/modules/auxil/psi_d_serial_mod.f90 similarity index 100% rename from base/modules/aux/psi_d_serial_mod.f90 rename to base/modules/auxil/psi_d_serial_mod.f90 diff --git a/base/modules/aux/psi_i_serial_mod.f90 b/base/modules/auxil/psi_i_serial_mod.f90 similarity index 100% rename from base/modules/aux/psi_i_serial_mod.f90 rename to base/modules/auxil/psi_i_serial_mod.f90 diff --git a/base/modules/aux/psi_s_serial_mod.f90 b/base/modules/auxil/psi_s_serial_mod.f90 similarity index 100% rename from base/modules/aux/psi_s_serial_mod.f90 rename to base/modules/auxil/psi_s_serial_mod.f90 diff --git a/base/modules/aux/psi_serial_mod.f90 b/base/modules/auxil/psi_serial_mod.f90 similarity index 100% rename from base/modules/aux/psi_serial_mod.f90 rename to base/modules/auxil/psi_serial_mod.f90 diff --git a/base/modules/aux/psi_z_serial_mod.f90 b/base/modules/auxil/psi_z_serial_mod.f90 similarity index 100% rename from base/modules/aux/psi_z_serial_mod.f90 rename to base/modules/auxil/psi_z_serial_mod.f90 diff --git a/base/modules/cutil.c b/base/modules/cutil.c index 13b7e1e94..11b815051 100644 --- a/base/modules/cutil.c +++ b/base/modules/cutil.c @@ -1,4 +1,6 @@ +#if ! (defined(_WIN32) || defined(WIN32)) #include +#endif #include #include #include "psb_internals.h" @@ -15,7 +17,7 @@ #ifdef UpperUnderscore #define psi_c_diffadd PSI_C_DIFFADD_ #endif -#ifdef UpperDoubleUnderscore +#ifdef UpperDoubleUnderscore #define psi_c_diffadd PSI_C_DIFFADD__ #endif #ifdef UpperCase diff --git a/base/modules/fakempi.c b/base/modules/fakempi.c index 9bf817334..c37afe271 100644 --- a/base/modules/fakempi.c +++ b/base/modules/fakempi.c @@ -1,4 +1,10 @@ +#if (defined(_WIN32) || defined(WIN32)) +#include +#include +#include +#else #include +#endif #include #include #include "psb_internals.h" @@ -102,8 +108,17 @@ #define mpi_complex 5 #define mpi_double_complex 6 -double mpi_wtime() +double mpi_wtime() { +#if defined(WIN32) || defined(_WIN32) + LARGE_INTEGER tim, freq; + double seconds; + + QueryPerformanceCounter(&tim); + QueryPerformanceFrequency(&freq); + seconds = (double)(tim.QuadPart) / (double)(freq.QuadPart); + return(seconds); +#else struct timeval tt; struct timezone tz; double temp; @@ -114,6 +129,7 @@ double mpi_wtime() temp = ((double)tt.tv_sec) + ((double)tt.tv_usec)*1.0e-6; } return(temp); +#endif } diff --git a/base/modules/psi_comm_buffers_mod.F90 b/base/modules/psi_comm_buffers_mod.F90 index fb418bb41..7c761226a 100644 --- a/base/modules/psi_comm_buffers_mod.F90 +++ b/base/modules/psi_comm_buffers_mod.F90 @@ -1,9 +1,9 @@ -! +! ! Parallel Sparse BLAS version 3.5 ! (C) Copyright 2006, 2010, 2015, 2017 ! Salvatore Filippone Cranfield University ! Alfredo Buttari CNRS-IRIT, Toulouse -! +! ! Redistribution and use in source and binary forms, with or without ! modification, are permitted provided that the following conditions ! are met: @@ -15,7 +15,7 @@ ! 3. The name of the PSBLAS group or the names of its contributors may ! not be used to endorse or promote products derived from this ! software without specific 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 @@ -27,9 +27,9 @@ ! 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. -! -! -#if defined(SERIAL_MPI) +! +! +#ifdef SERIAL_MPI ! Provide a fake mpi module just to keep the compiler(s) happy. module mpi use psb_const_mod @@ -40,17 +40,17 @@ module mpi integer(psb_mpik_), parameter :: mpi_integer8 = 2 integer(psb_mpik_), parameter :: mpi_real = 3 integer(psb_mpik_), parameter :: mpi_double_precision = 4 - integer(psb_mpik_), parameter :: mpi_complex = 5 - integer(psb_mpik_), parameter :: mpi_double_complex = 6 + integer(psb_mpik_), parameter :: mpi_complex = 5 + integer(psb_mpik_), parameter :: mpi_double_complex = 6 integer(psb_mpik_), parameter :: mpi_character = 7 integer(psb_mpik_), parameter :: mpi_logical = 8 integer(psb_mpik_), parameter :: mpi_integer2 = 9 integer(psb_mpik_), parameter :: mpi_comm_null = -1 integer(psb_mpik_), parameter :: mpi_comm_world = 1 - + real(psb_dpk_), external :: mpi_wtime end module mpi -#endif +#endif module psi_comm_buffers_mod use psb_const_mod @@ -69,7 +69,7 @@ module psi_comm_buffers_mod type psb_buffer_node integer(psb_mpik_) :: request - integer(psb_mpik_) :: icontxt + integer(psb_mpik_) :: icontxt integer(psb_mpik_) :: buffer_type integer(psb_ipk_), allocatable :: intbuf(:) integer(psb_long_int_k_), allocatable :: int8buf(:) @@ -112,25 +112,25 @@ module psi_comm_buffers_mod module procedure psi_i2snd end interface #endif - + contains subroutine psb_init_queue(mesg_queue,info) - implicit none + implicit none type(psb_buffer_queue), intent(inout) :: mesg_queue integer(psb_ipk_), intent(out) :: info info = 0 if ((.not.associated(mesg_queue%head)).and.& - & (.not.associated(mesg_queue%tail))) then + & (.not.associated(mesg_queue%tail))) then ! Nothing to do return end if if ((.not.associated(mesg_queue%head)).or.& - & (.not.associated(mesg_queue%tail))) then + & (.not.associated(mesg_queue%tail))) then ! If we are here one is associated, the other is not. - ! This is impossible. + ! This is impossible. info = -1 write(psb_err_unit,*) 'Wrong status on init ' return @@ -142,12 +142,12 @@ subroutine psb_wait_buffer(node, info) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif type(psb_buffer_node), intent(inout) :: node - integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), intent(out) :: info integer(psb_mpik_) :: status(mpi_status_size),minfo minfo = mpi_success call mpi_wait(node%request,status,minfo) @@ -158,13 +158,13 @@ subroutine psb_test_buffer(node, flag, info) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif type(psb_buffer_node), intent(inout) :: node logical, intent(out) :: flag - integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), intent(out) :: info integer(psb_mpik_) :: status(mpi_status_size), minfo minfo = mpi_success #if defined(SERIAL_MPI) @@ -174,7 +174,7 @@ subroutine psb_test_buffer(node, flag, info) #endif info=minfo end subroutine psb_test_buffer - + subroutine psb_close_context(mesg_queue,icontxt) type(psb_buffer_queue), intent(inout) :: mesg_queue @@ -183,10 +183,10 @@ subroutine psb_close_context(mesg_queue,icontxt) type(psb_buffer_node), pointer :: node, nextnode node => mesg_queue%head - do + do if (.not.associated(node)) exit nextnode => node%next - if (node%icontxt == icontxt) then + if (node%icontxt == icontxt) then call psb_wait_buffer(node,info) call psb_delete_node(mesg_queue,node) end if @@ -198,9 +198,9 @@ subroutine psb_close_all_context(mesg_queue) type(psb_buffer_queue), intent(inout) :: mesg_queue type(psb_buffer_node), pointer :: node, nextnode integer(psb_ipk_) :: info - + node => mesg_queue%head - do + do if (.not.associated(node)) exit nextnode => node%next call psb_wait_buffer(node,info) @@ -214,8 +214,8 @@ subroutine psb_delete_node(mesg_queue,node) type(psb_buffer_queue), intent(inout) :: mesg_queue type(psb_buffer_node), pointer :: node type(psb_buffer_node), pointer :: prevnode - - if (.not.associated(node)) then + + if (.not.associated(node)) then return end if prevnode => node%prev @@ -224,7 +224,7 @@ subroutine psb_delete_node(mesg_queue,node) if (associated(prevnode)) prevnode%next => node%next if (associated(node%next)) node%next%prev => prevnode deallocate(node) - + end subroutine psb_delete_node subroutine psb_insert_node(mesg_queue,node) @@ -234,7 +234,7 @@ subroutine psb_insert_node(mesg_queue,node) node%next => null() node%prev => null() if ((.not.associated(mesg_queue%head)).and.& - & (.not.associated(mesg_queue%tail))) then + & (.not.associated(mesg_queue%tail))) then mesg_Queue%head => node mesg_queue%tail => node return @@ -250,13 +250,13 @@ subroutine psb_test_nodes(mesg_queue) type(psb_buffer_node), pointer :: node, nextnode integer(psb_ipk_) :: info logical :: flag - + node => mesg_queue%head - do + do if (.not.associated(node)) exit nextnode => node%next call psb_test_buffer(node,flag,info) - if (flag) then + if (flag) then call psb_delete_node(mesg_queue,node) end if node => nextnode @@ -270,14 +270,14 @@ end subroutine psb_test_nodes ! to a node in the mesg queue, then it is sent. ! Thus the calling process should guarantee that ! the buffer is dispensable, i.e. the user data - ! has already been copied. + ! has already been copied. ! ! !!!!!!!!!!!!!!!!! subroutine psi_isnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -287,16 +287,16 @@ subroutine psi_isnd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_ipk_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_int_type call move_alloc(buffer,node%intbuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -304,7 +304,7 @@ subroutine psi_isnd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) end subroutine psi_isnd @@ -314,7 +314,7 @@ subroutine psi_i4snd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -324,24 +324,24 @@ subroutine psi_i4snd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_mpik_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_int4_type call move_alloc(buffer,node%int4buf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if call mpi_isend(node%int4buf,size(node%int4buf),psb_mpi_def_integer,& & dest,tag,icontxt,node%request,minfo) - info = minfo + info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) end subroutine psi_i4snd @@ -352,7 +352,7 @@ subroutine psi_i8snd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -362,24 +362,24 @@ subroutine psi_i8snd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_ipk_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_int8_type call move_alloc(buffer,node%int8buf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if call mpi_isend(node%int8buf,size(node%int8buf),psb_mpi_lng_integer,& & dest,tag,icontxt,node%request,minfo) - info = minfo + info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) end subroutine psi_i8snd @@ -390,7 +390,7 @@ subroutine psi_i2snd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -400,16 +400,16 @@ subroutine psi_i2snd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_ipk_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_int2_type call move_alloc(buffer,node%int2buf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -417,7 +417,7 @@ subroutine psi_i2snd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) end subroutine psi_i2snd @@ -427,7 +427,7 @@ subroutine psi_ssnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -439,14 +439,14 @@ subroutine psi_ssnd(icontxt,tag,dest,buffer,mesg_queue) integer(psb_mpik_) :: minfo allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_real_type call move_alloc(buffer,node%realbuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -454,16 +454,16 @@ subroutine psi_ssnd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) - + end subroutine psi_ssnd subroutine psi_dsnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -475,14 +475,14 @@ subroutine psi_dsnd(icontxt,tag,dest,buffer,mesg_queue) integer(psb_mpik_) :: minfo allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_double_type call move_alloc(buffer,node%doublebuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -490,16 +490,16 @@ subroutine psi_dsnd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) - + end subroutine psi_dsnd - + subroutine psi_csnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -511,31 +511,31 @@ subroutine psi_csnd(icontxt,tag,dest,buffer,mesg_queue) integer(psb_mpik_) :: minfo allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_complex_type call move_alloc(buffer,node%complexbuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if call mpi_isend(node%complexbuf,size(node%complexbuf),psb_mpi_c_spk_,& & dest,tag,icontxt,node%request,minfo) - info = minfo + info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) - + end subroutine psi_csnd subroutine psi_zsnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -545,16 +545,16 @@ subroutine psi_zsnd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_ipk_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_dcomplex_type call move_alloc(buffer,node%dcomplbuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -562,9 +562,9 @@ subroutine psi_zsnd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) - + end subroutine psi_zsnd @@ -572,7 +572,7 @@ subroutine psi_lsnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -582,16 +582,16 @@ subroutine psi_lsnd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_ipk_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_logical_type call move_alloc(buffer,node%logbuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -599,9 +599,9 @@ subroutine psi_lsnd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) - + end subroutine psi_lsnd @@ -609,7 +609,7 @@ subroutine psi_hsnd(icontxt,tag,dest,buffer,mesg_queue) #ifdef MPI_MOD use mpi #endif - implicit none + implicit none #ifdef MPI_H include 'mpif.h' #endif @@ -619,16 +619,16 @@ subroutine psi_hsnd(icontxt,tag,dest,buffer,mesg_queue) type(psb_buffer_node), pointer :: node integer(psb_ipk_) :: info integer(psb_mpik_) :: minfo - + allocate(node, stat=info) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if node%icontxt = icontxt node%buffer_type = psb_char_type call move_alloc(buffer,node%charbuf) - if (info /= 0) then + if (info /= 0) then write(psb_err_unit,*) 'Fatal memory error inside communication subsystem' return end if @@ -636,9 +636,9 @@ subroutine psi_hsnd(icontxt,tag,dest,buffer,mesg_queue) & dest,tag,icontxt,node%request,minfo) info = minfo call psb_insert_node(mesg_queue,node) - + call psb_test_nodes(mesg_queue) - + end subroutine psi_hsnd diff --git a/base/psblas/psb_sxdot.f90 b/base/psblas/psb_sxdot.f90 deleted file mode 100644 index aafe474f9..000000000 --- a/base/psblas/psb_sxdot.f90 +++ /dev/null @@ -1,589 +0,0 @@ -! -! Parallel Sparse BLAS version 3.5 -! (C) Copyright 2006, 2010, 2015, 2017 -! Salvatore Filippone Cranfield University -! Alfredo Buttari CNRS-IRIT, Toulouse -! -! 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 PSBLAS group or the names of its contributors may -! not be used to endorse or promote products derived from this -! software without specific 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 PSBLAS 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: psb_sdot.f90 -! -! Function: psb_sdot -! psb_sdot forms the dot product of two distributed vectors, -! -! dot := sub( X )**T * sub( Y ) -! -! where sub( X ) denotes X(:,JX) -! -! sub( Y ) denotes Y(:,JY). -! -! Arguments: -! x(:,:) - real The input vector containing the entries of ( X ). -! y(:,:) - real The input vector containing the entries of ( Y ). -! desc_a - type(psb_desc_type). The communication descriptor. -! info - integer. Return code -! jx - integer(optional). The column offset for sub( X ). -! jy - integer(optional). The column offset for sub( Y ). -! -function psb_sxdot(x, y,desc_a, info, jx, jy) - use psb_desc_mod - use psb_check_mod - use psb_error_mod - use psb_penv_mod - implicit none - - real(psb_spk_), intent(in) :: x(:,:), y(:,:) - type(psb_desc_type), intent(in) :: desc_a - integer(psb_ipk_), intent(in), optional :: jx, jy - integer(psb_ipk_), intent(out) :: info - real(psb_dpk_) :: psb_sxdot - - ! locals - integer(psb_ipk_) :: ictxt, np, me, idx, ndm,& - & err_act, iix, jjx, ix, ijx, iy, ijy, iiy, jjy, i, m - real(psb_dpk_) :: dot_local - real(psb_dpk_) :: sxdot - character(len=20) :: name, ch_err - - name='psb_sdot' - if(psb_get_errstatus() /= 0) return - info=psb_success_ - call psb_erractionsave(err_act) - - ictxt=desc_a%get_context() - call psb_info(ictxt, me, np) - if (np == -ione) then - info = psb_err_context_error_ - call psb_errpush(info,name) - goto 9999 - endif - - ix = ione - if (present(jx)) then - ijx = jx - else - ijx = ione - endif - - iy = ione - if (present(jy)) then - ijy = jy - else - ijy = ione - endif - - if(ijx /= ijy) then - info=3050 - call psb_errpush(info,name) - goto 9999 - end if - - m = desc_a%get_global_rows() - - ! check vector correctness - call psb_chkvect(m,ione,size(x,1),ix,ijx,desc_a,info,iix,jjx) - if (info == psb_success_) & - & call psb_chkvect(m,ione,size(y,1),iy,ijy,desc_a,info,iiy,jjy) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - ch_err='psb_chkvect' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - - if ((iix /= ione).or.(iiy /= ione)) then - info=psb_err_ix_n1_iy_n1_unsupported_ - call psb_errpush(info,name) - goto 9999 - end if - - if(m /= 0) then - if(desc_a%get_local_rows() > 0) then - dot_local = sxdot(desc_a%get_local_rows(),& - & x(iix,jjx),ione,y(iiy,jjy),ione) - ! adjust dot_local because overlapped elements are computed more than once - do i=1,size(desc_a%ovrlap_elem,1) - idx = desc_a%ovrlap_elem(i,1) - ndm = desc_a%ovrlap_elem(i,2) - dot_local = dot_local - (real(ndm-1)/real(ndm))*(x(idx,jjx)*y(idx,jjy)) - end do - else - dot_local=0.0 - end if - else - dot_local=0.0 - end if - - ! compute global sum - call psb_sum(ictxt, dot_local) - - psb_sxdot = dot_local - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - - if (err_act == psb_act_abort_) then - call psb_error(ictxt) - return - end if - return -end function psb_sxdot - - - - -!!$ -!!$ Parallel Sparse BLAS version 3.5 -!!$ (C) Copyright 2006, 2010, 2015, 2017 -!!$ Salvatore Filippone University of Rome Tor Vergata -!!$ Alfredo Buttari CNRS-IRIT, Toulouse -!!$ -!!$ 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 PSBLAS group or the names of its contributors may -!!$ not be used to endorse or promote products derived from this -!!$ software without specific 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 PSBLAS 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. -!!$ -!!$ -! -! Function: psb_sdotv -! psb_sdotv forms the dot product of two distributed vectors, -! -! dot := X**T * Y -! -! Arguments: -! x(:) - real The input vector containing the entries of X. -! y(:) - real The input vector containing the entries of Y. -! desc_a - type(psb_desc_type). The communication descriptor. -! info - integer. Return code -! -function psb_sxdotv(x, y,desc_a, info) - use psb_desc_mod - use psb_check_mod - use psb_error_mod - use psb_penv_mod - implicit none - - real(psb_spk_), intent(in) :: x(:), y(:) - type(psb_desc_type), intent(in) :: desc_a - integer(psb_ipk_), intent(out) :: info - real(psb_dpk_) :: psb_sxdotv - - ! locals - integer(psb_ipk_) :: ictxt, np, me, idx, ndm,& - & err_act, iix, jjx, ix, jx, iy, jy, iiy, jjy, i, m - real(psb_dpk_) :: dot_local - real(psb_dpk_) :: sxdot - character(len=20) :: name, ch_err - - name='psb_sdot' - if(psb_get_errstatus() /= 0) return - info=psb_success_ - call psb_erractionsave(err_act) - - ictxt=desc_a%get_context() - - call psb_info(ictxt, me, np) - if (np == -ione) then - info = psb_err_context_error_ - call psb_errpush(info,name) - goto 9999 - endif - - ix = ione - iy = ione - jx = ione - jy = ione - m = desc_a%get_global_rows() - - ! check vector correctness - call psb_chkvect(m,ione,size(x,1),ix,jx,desc_a,info,iix,jjx) - if (info == psb_success_) & - & call psb_chkvect(m,ione,size(y,1),iy,jy,desc_a,info,iiy,jjy) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - ch_err='psb_chkvect' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - - if ((iix /= ione).or.(iiy /= ione)) then - info=psb_err_ix_n1_iy_n1_unsupported_ - call psb_errpush(info,name) - goto 9999 - end if - - if(m /= 0) then - if(desc_a%get_local_rows() > 0) then - dot_local = sxdot(desc_a%get_local_rows(),& - & x,ione,y,ione) - ! adjust dot_local because overlapped elements are computed more than once - do i=1,size(desc_a%ovrlap_elem,1) - idx = desc_a%ovrlap_elem(i,1) - ndm = desc_a%ovrlap_elem(i,2) - dot_local = dot_local - (real(ndm-1)/real(ndm))*(x(idx)*y(idx)) - end do - else - dot_local=0.0 - end if - else - dot_local=0.0 - end if - - ! compute global sum - call psb_sum(ictxt, dot_local) - - psb_sxdotv = dot_local - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - - if (err_act == psb_act_abort_) then - call psb_error(ictxt) - return - end if - return -end function psb_sxdotv - -function sxdot(n,x,ix,y,iy) - use psb_const_mod - real(psb_dpk_) :: sxdot - integer(psb_ipk_) :: n,ix,iy - real(psb_spk_) :: x(*),y(*) - real(psb_dpk_) :: tmp - integer(psb_ipk_) :: i - - if ((ix /= 1).or.(iy /= 1)) then - write(psb_err_unit,*) 'WARNING unimplemented case in SXDOT' - sxdot=dzero - return - end if - - tmp = dzero - do i=1, n - tmp = tmp + dble(x(i))*dble(y(i)) - end do - sxdot = tmp - -end function sxdot - - -!!$ -!!$ Parallel Sparse BLAS version 3.5 -!!$ (C) Copyright 2006, 2010, 2015, 2017 -!!$ Salvatore Filippone University of Rome Tor Vergata -!!$ Alfredo Buttari CNRS-IRIT, Toulouse -!!$ -!!$ 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 PSBLAS group or the names of its contributors may -!!$ not be used to endorse or promote products derived from this -!!$ software without specific 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 PSBLAS 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. -!!$ -!!$ -! -! Subroutine: psb_sdotvs -! psb_sdotvs forms the dot product of two distributed vectors, -! -! dot := X**T * Y -! -! Arguments: -! res - real The result. -! x(:) - real The input vector containing the entries of X. -! y(:) - real The input vector containing the entries of Y. -! desc_a - type(psb_desc_type). The communication descriptor. -! info - integer. Return code -! -subroutine psb_sxdotvs(res, x, y,desc_a, info) - use psb_base_mod, psb_protect_name => psb_sxdotvs - implicit none - - real(psb_spk_), intent(in) :: x(:), y(:) - real(psb_dpk_), intent(out) :: res - type(psb_desc_type), intent(in) :: desc_a - integer(psb_ipk_), intent(out) :: info - - ! locals - integer(psb_ipk_) :: ictxt, np, me, idx, ndm,& - & err_act, iix, jjx, ix, iy, iiy, jjy, i, m - real(psb_dpk_) :: dot_local - real(psb_dpk_) :: sxdot - character(len=20) :: name, ch_err - - name='psb_sdot' - if(psb_get_errstatus() /= 0) return - info=psb_success_ - call psb_erractionsave(err_act) - - ictxt=desc_a%get_context() - - call psb_info(ictxt, me, np) - if (np == -ione) then - info = psb_err_context_error_ - call psb_errpush(info,name) - goto 9999 - endif - - ix = ione - iy = ione - m = desc_a%get_global_rows() - ! check vector correctness - call psb_chkvect(m,ione,size(x,1),ix,ix,desc_a,info,iix,jjx) - if (info == psb_success_) & - & call psb_chkvect(m,ione,size(y,1),iy,iy,desc_a,info,iiy,jjy) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - ch_err='psb_chkvect' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - - if ((iix /= ione).or.(iiy /= ione)) then - info=psb_err_ix_n1_iy_n1_unsupported_ - call psb_errpush(info,name) - goto 9999 - end if - - if(m /= 0) then - if(desc_a%get_local_rows() > 0) then - dot_local = sxdot(desc_a%get_local_rows(),& - & x,ione,y,ione) - ! adjust dot_local because overlapped elements are computed more than once - do i=1,size(desc_a%ovrlap_elem,1) - idx = desc_a%ovrlap_elem(i,1) - ndm = desc_a%ovrlap_elem(i,2) - dot_local = dot_local - (real(ndm-1)/real(ndm))*(x(idx)*y(idx)) - end do - else - dot_local=0.0 - end if - else - dot_local=0.0 - end if - - ! compute global sum - call psb_sum(ictxt, dot_local) - - res = dot_local - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - - if (err_act == psb_act_abort_) then - call psb_error(ictxt) - return - end if - return -end subroutine psb_sxdotvs - - - - -!!$ -!!$ Parallel Sparse BLAS version 3.5 -!!$ (C) Copyright 2006, 2010, 2015, 2017 -!!$ Salvatore Filippone University of Rome Tor Vergata -!!$ Alfredo Buttari CNRS-IRIT, Toulouse -!!$ -!!$ 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 PSBLAS group or the names of its contributors may -!!$ not be used to endorse or promote products derived from this -!!$ software without specific 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 PSBLAS 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. -!!$ -!!$ -! -! Subroutine: psb_smdots -! psb_smdots forms the dot product of multiple distributed vectors, -! -! res(i) := ( X(:,i) )**T * ( Y(:,i) ) -! -! Arguments: -! res(:) - real. The result. -! x(:,:) - real The input vector containing the entries of ( X ). -! y(:,:) - real The input vector containing the entries of ( Y ). -! desc_a - type(psb_desc_type). The communication descriptor. -! info - integer. Return code -! -subroutine psb_sxmdots(res, x, y, desc_a, info) - use psb_base_mod, psb_protect_name => psb_sxmdots - implicit none - - real(psb_spk_), intent(in) :: x(:,:), y(:,:) - real(psb_dpk_), intent(out) :: res(:) - type(psb_desc_type), intent(in) :: desc_a - integer(psb_ipk_), intent(out) :: info - - ! locals - integer(psb_ipk_) :: ictxt, np, me, idx, ndm,& - & err_act, iix, jjx, ix, iy, iiy, jjy, i, m, j, k - real(psb_dpk_),allocatable :: dot_local(:) - real(psb_dpk_) :: sxdot - character(len=20) :: name, ch_err - - name='psb_smdots' - if(psb_get_errstatus() /= 0) return - info=psb_success_ - call psb_erractionsave(err_act) - - ictxt=desc_a%get_context() - - call psb_info(ictxt, me, np) - if (np == -ione) then - info = psb_err_context_error_ - call psb_errpush(info,name) - goto 9999 - endif - - ix = ione - iy = ione - - m = desc_a%get_global_rows() - - ! check vector correctness - call psb_chkvect(m,ione,size(x,1),ix,ix,desc_a,info,iix,jjx) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - ch_err='psb_chkvect' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - call psb_chkvect(m,ione,size(y,1),iy,iy,desc_a,info,iiy,jjy) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - ch_err='psb_chkvect' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - - if ((ix /= ione).or.(iy /= ione)) then - info=psb_err_ix_n1_iy_n1_unsupported_ - call psb_errpush(info,name) - goto 9999 - end if - - k = min(size(x,2),size(y,2)) - allocate(dot_local(k)) - - if(m /= 0) then - if(desc_a%get_local_rows() > 0) then - do j=1,k - dot_local(j) = sxdot(desc_a%get_local_rows(),& - & x(1,j),ione,y(1,j),ione) - ! adjust dot_local because overlapped elements are computed more than once - end do - do i=1,size(desc_a%ovrlap_elem,1) - idx = desc_a%ovrlap_elem(i,1) - ndm = desc_a%ovrlap_elem(i,2) - dot_local(1:k) = dot_local(1:k) - (real(ndm-1)/real(ndm))*(x(idx,1:k)*y(idx,1:k)) - end do - else - dot_local(:)=0.0 - end if - else - dot_local(:)=0.0 - end if - - ! compute global sum - call psb_sum(ictxt, dot_local(1:k)) - - res(1:k) = dot_local(1:k) - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - - if (err_act == psb_act_abort_) then - call psb_error(ictxt) - return - end if - return -end subroutine psb_sxmdots diff --git a/cmake/CapitalizeString.cmake b/cmake/CapitalizeString.cmake new file mode 100644 index 000000000..ae50b4745 --- /dev/null +++ b/cmake/CapitalizeString.cmake @@ -0,0 +1,7 @@ +function(CapitalizeString string output_variable) + string(TOUPPER "${string}" _upper_string) + string(TOLOWER "${string}" _lower_string) + string(SUBSTRING "${_upper_string}" 0 1 _start) + string(SUBSTRING "${_lower_string}" 1 -1 _end) + set(${output_variable} "${_start}${_end}" PARENT_SCOPE) +endfunction() diff --git a/cmake/CheckOutOfSourceBuild.cmake b/cmake/CheckOutOfSourceBuild.cmake new file mode 100644 index 000000000..ad1769522 --- /dev/null +++ b/cmake/CheckOutOfSourceBuild.cmake @@ -0,0 +1,21 @@ +#-------------------------- +# Prohibit in-source builds +#-------------------------- +if ("${CMAKE_CURRENT_SOURCE_DIR}" STREQUAL "${CMAKE_CURRENT_BINARY_DIR}") + message(FATAL_ERROR "ERROR! " + "CMAKE_CURRENT_SOURCE_DIR=${CMAKE_CURRENT_SOURCE_DIR}" + " == CMAKE_CURRENT_BINARY_DIR=${CMAKE_CURRENT_BINARY_DIR}" + "\nThis archive does not support in-source builds:\n" + "You must now delete the CMakeCache.txt file and the CMakeFiles/ directory under " + "the 'src' source directory or you will not be able to configure correctly!" + "\nYou must now run something like:\n" + " $ rm -r CMakeCache.txt CMakeFiles/" + "\n" + "Please create a directory outside the ${CMAKE_PROJECT_NAME} source tree and build under that outside directory " + "in a manner such as\n" + " $ mkdir build\n" + " $ cd build\n" + " $ CC=gcc FC=gfortran cmake -DBUILD_TYPE=Release -DCMAKE_INSTALL_PREFIX=/path/to/install/dir /path/to/psblas3/src/dir \n" + "\nsubstituting the appropriate syntax for your shell (the above line assumes the bash shell)." + ) +endif() diff --git a/cmake/FindMETIS.cmake b/cmake/FindMETIS.cmake new file mode 100644 index 000000000..21f95f572 --- /dev/null +++ b/cmake/FindMETIS.cmake @@ -0,0 +1,95 @@ +if (METIS_INCLUDES AND METIS_LIBRARIES) + set(METIS_FIND_QUIETLY TRUE) +endif (METIS_INCLUDES AND METIS_LIBRARIES) + +if( DEFINED ENV{METISDIR} ) + if( NOT DEFINED METIS_ROOT ) + set(METIS_ROOT "$ENV{METISDIR}") + endif() +endif() + +if( (DEFINED ENV{METIS_ROOT}) OR (DEFINED METIS_ROOT) ) + if( NOT DEFINED METIS_ROOT) + set(METIS_ROOT "$ENV{METIS_ROOT}") + endif() + set(METIS_HINTS "${METIS_ROOT}") +endif() + +find_path(METIS_INCLUDES + NAMES + metis.h + HINTS + ${METIS_HINTS} + PATHS + "${INCLUDE_INSTALL_DIR}" + /usr/local/opt + /usr/local + PATH_SUFFIXES + include + ) + +if(METIS_INCLUDES) + foreach(include IN_LISTS METIS_INCLUDES) + get_filename_component(mts_include_dir "${include}" DIRECTORY) + get_filename_component(mts_abs_include_dir "${mts_include_dir}" ABSOLUTE) + get_filename_component(new_mts_hint "${include_dir}/.." ABSOLUTE ) + list(APPEND METIS_HINTS "${new_mts_hint}") + break() + endforeach() +endif() + +if(METIS_HINTS) + list(REMOVE_DUPLICATES METIS_HINTS) +endif() + +macro(_metis_check_version) + file(READ "${METIS_INCLUDES}/metis.h" _metis_version_header) + + string(REGEX MATCH "define[ \t]+METIS_VER_MAJOR[ \t]+([0-9]+)" _metis_major_version_match "${_metis_version_header}") + set(METIS_MAJOR_VERSION "${CMAKE_MATCH_1}") + string(REGEX MATCH "define[ \t]+METIS_VER_MINOR[ \t]+([0-9]+)" _metis_minor_version_match "${_metis_version_header}") + set(METIS_MINOR_VERSION "${CMAKE_MATCH_1}") + string(REGEX MATCH "define[ \t]+METIS_VER_SUBMINOR[ \t]+([0-9]+)" _metis_subminor_version_match "${_metis_version_header}") + set(METIS_SUBMINOR_VERSION "${CMAKE_MATCH_1}") + if(NOT METIS_MAJOR_VERSION) + message(STATUS "Could not determine Metis version. Assuming version 4.0.0") + set(METIS_VERSION 4.0.0) + else() + set(METIS_VERSION ${METIS_MAJOR_VERSION}.${METIS_MINOR_VERSION}.${METIS_SUBMINOR_VERSION}) + endif() + if(${METIS_VERSION} VERSION_LESS ${Metis_FIND_VERSION}) + set(METIS_VERSION_OK FALSE) + else() + set(METIS_VERSION_OK TRUE) + endif() + + if(NOT METIS_VERSION_OK) + message(STATUS "Metis version ${METIS_VERSION} found in ${METIS_INCLUDES}, " + "but at least version ${Metis_FIND_VERSION} is required") + endif(NOT METIS_VERSION_OK) +endmacro(_metis_check_version) + +if(METIS_INCLUDES AND Metis_FIND_VERSION) + _metis_check_version() +else() + set(METIS_VERSION_OK TRUE) +endif() + + +find_library(METIS_LIBRARIES metis + HINTS + ${METIS_HINTS} + PATHS + "${LIB_INSTALL_DIR}" + /usr/local/ + /usr/local/opt + PATH_SUFFIXES + lib + lib64 + metis/lib) + +include(FindPackageHandleStandardArgs) +find_package_handle_standard_args(METIS DEFAULT_MSG + METIS_INCLUDES METIS_LIBRARIES METIS_VERSION_OK) + +mark_as_advanced(METIS_INCLUDES METIS_LIBRARIES) diff --git a/cmake/makeDist.cmake b/cmake/makeDist.cmake new file mode 100644 index 000000000..f0b465660 --- /dev/null +++ b/cmake/makeDist.cmake @@ -0,0 +1,79 @@ +# CMake file to be called in script mode (${CMAKE_COMMAND} -P ) to +# Generate a source archive release asset from add_custom_command +# +# See SourceDistTarget.cmake + +if(NOT CMAKE_ARGV3) + message(FATAL_ERROR "Must pass the top level src dir to ${CMAKE_ARGV2} as the first argument") +endif() + +if(NOT CMAKE_ARGV4) + message(FATAL_ERROR "Must pass the top level src dir to ${CMAKE_ARGV2} as the second argument") +endif() + +find_package(Git) +if(NOT GIT_FOUND) + message( FATAL_ERROR "You can't create a source archive release asset without git!") +endif() + +execute_process(COMMAND "${GIT_EXECUTABLE}" describe --always + RESULT_VARIABLE git_status + OUTPUT_VARIABLE git_version + WORKING_DIRECTORY "${CMAKE_ARGV3}" + OUTPUT_STRIP_TRAILING_WHITESPACE) +if(NOT (git_status STREQUAL "0")) + message( FATAL_ERROR "git describe --always failed with exit status: ${git_status} and message: +${git_version}") +endif() + +set(archive "PSBLAS-${git_version}") +set(l_archive "PSBLAS-${git_version}") +set(release_asset "${CMAKE_ARGV4}/${archive}.tar.gz") +execute_process( + COMMAND "${GIT_EXECUTABLE}" archive "--prefix=${archive}/" -o "${release_asset}" "${git_version}" + RESULT_VARIABLE git_status + OUTPUT_VARIABLE git_output + WORKING_DIRECTORY "${CMAKE_ARGV3}" + OUTPUT_STRIP_TRAILING_WHITESPACE) + +if(NOT (git_status STREQUAL "0")) + message( FATAL_ERROR "git archive ... failed with exit status: ${git_status} and message: +${git_output}") +else() + message( STATUS "Source code release asset created from `git archive`: ${release_asset}") +endif() + +file(SHA256 "${release_asset}" tarball_sha256) +set(sha256_checksum "${tarball_sha256} ${archive}.tar.gz") +configure_file("${CMAKE_ARGV3}/cmake/PSBLAS-VER-SHA256.txt.in" + "${CMAKE_ARGV4}/${l_archive}-SHA256.txt" + @ONLY) +message( STATUS + "SHA 256 checksum of release tarball written out as: ${CMAKE_ARGV4}/${l_archive}-SHA256.txt" ) + +find_program(GPG_EXECUTABLE + gpg + DOC "Location of GnuPG (gpg) executable") + +if(GPG_EXECUTABLE) + execute_process( + COMMAND "${GPG_EXECUTABLE}" --armor --detach-sign --comment "@gpg_comment@" "${CMAKE_ARGV4}/${l_archive}-SHA256.txt" + RESULT_VARIABLE gpg_status + OUTPUT_VARIABLE gpg_output + WORKING_DIRECTORY "${CMAKE_ARGV4}") + if(NOT (gpg_status STREQUAL "0")) + message( WARNING "GPG signing of ${CMAKE_ARGV4}/${l_archive}-SHA256.txt appears to have failed +with status: ${gpg_status} and output: ${gpg_output}") + else() + configure_file("${CMAKE_ARGV3}/cmake/PSBLAS-VER-SHA256.txt.asc.in" + "${CMAKE_ARGV4}/${l_archive}-GPG.comment" + @ONLY) + file(READ "${CMAKE_ARGV4}/${l_archive}-GPG.comment" gpg_comment) + configure_file("${CMAKE_ARGV4}/${l_archive}-SHA256.txt.asc" + "${CMAKE_ARGV4}/${l_archive}-SHA256.txt.asc.out" + @ONLY) + file(RENAME "${CMAKE_ARGV4}/${l_archive}-SHA256.txt.asc.out" + "${CMAKE_ARGV4}/${l_archive}-SHA256.txt.asc") + message(STATUS "GPG signed SHA256 checksum created: ${CMAKE_ARGV4}/${l_archive}-SHA256.txt.asc") + endif() +endif() diff --git a/cmake/pkg/psblasConfig.cmake.in b/cmake/pkg/psblasConfig.cmake.in new file mode 100644 index 000000000..3fe60345b --- /dev/null +++ b/cmake/pkg/psblasConfig.cmake.in @@ -0,0 +1,16 @@ +# Config file for the INSTALLED package +# Allow other CMake projects to find this package if it is installed +# Requires the use of the standard CMake module CMakePackageConfigHelpers + +set ( @CMAKE_PROJECT_NAME@_VERSION @VERSION@ ) + +###@COMPILER_CONSISTENCY_CHECK@ + +@PACKAGE_INIT@ + +# Provide the targets +set_and_check ( @PACKAGE_NAME@_CONFIG_INSTALL_DIR "@PACKAGE_EXPORT_INSTALL_DIR@" ) +include ( "${@PACKAGE_NAME@_CONFIG_INSTALL_DIR}/@PACKAGE_NAME@-targets.cmake" ) + +# Make the module files available via include +set_and_check ( @CMAKE_PROJECT_NAME@_INCLUDE_DIRS "@PACKAGE_INSTALL_MOD_DIR@" ) diff --git a/cmake/psblas-VER-SHA256.txt.asc.in b/cmake/psblas-VER-SHA256.txt.asc.in new file mode 100644 index 000000000..22cd9211b --- /dev/null +++ b/cmake/psblas-VER-SHA256.txt.asc.in @@ -0,0 +1,12 @@ +Mac users can use GPGTools - https://gpgtools.org +Comment: Download Izaak Beekman's GPG public key from your +Comment: trusted key server or from +Comment: https://izaakbeekman.com/izaak.pubkey.txt +Comment: Next add it to your GPG keyring, e.g., +Comment: `curl https://izaakbeekman.com/izaak.pubkey.txt | gpg --import` +Comment: Make sure you have verified that the release archive's +Comment: SHA256 checksum matches the provided +Comment: psblas-@git_version@-SHA256.txt and ensure that this file +Comment: and it's signature are in the same directory. Then +Comment: verify with: +Comment: `gpg --verify psblas-@git_version@-SHA256.txt.asc` diff --git a/cmake/psblas-VER-SHA256.txt.in b/cmake/psblas-VER-SHA256.txt.in new file mode 100644 index 000000000..1b3888c42 --- /dev/null +++ b/cmake/psblas-VER-SHA256.txt.in @@ -0,0 +1,3 @@ +# To verify cryptographic checksums `shasum -c psblas-@git_version@-SHA256.txt` on Mac OS X, or +# `sha256sum -c psblas-@git_version@-SHA256.txt` on Linux. +@sha256_checksum@ diff --git a/cmake/setVersion.cmake b/cmake/setVersion.cmake new file mode 100644 index 000000000..f79ee6195 --- /dev/null +++ b/cmake/setVersion.cmake @@ -0,0 +1,90 @@ +include(CMakeParseArguments) + +# Function to parse version info from git and/or .VERSION file +function(set_version) + set(options "") + set(oneValueArgs VERSION_VARIABLE GIT_DESCRIBE_VAR CUSTOM_VERSION_FILE CUSTOM_VERSION_REGEX ) + set(multiValueArgs "") + cmake_parse_arguments(set_version "${options}" "${oneValueArgs}" "${multiValueArgs}" ${ARGN}) + + # Algorithm: + # 1. Get first line of .VERSION file, which will be set via `git archive` so long as + # + # 2. If not a packaged release check if this is an active git repo + # 3. Get version info from `git describe` + # 4. First the most recent tag is fetched if available + # 5. Then the full `git describe` output is fetched + + + if(NOT set_version_CUSTOM_VERSION_REGEX) + set(_VERSION_REGEX "[vV]*[0-9]+\\.[0-9]+\\.[0-9]+") + else() + set(_VERSION_REGEX ${set_version_CUSTOM_VERSION_REGEX}) + endif() + if(NOT set_version_CUSTOM_VERSION_FILE) + set(_VERSION_FILE "${CMAKE_SOURCE_DIR}/.VERSION") + else() + set(_VERSION_FILE "${set_version_CUSTOM_VERSION_FILE}") + endif() + + file(STRINGS "${_VERSION_FILE}" first_line + LIMIT_COUNT 1 + ) + + string(REGEX MATCH ${_VERSION_REGEX} + _package_version "${first_line}") + + if((NOT (_package_version MATCHES ${_VERSION_REGEX})) AND (EXISTS "${CMAKE_SOURCE_DIR}/.git")) + message( STATUS "Build from git repository detected") + find_package(Git) + if(GIT_FOUND) + set(GIT_FOUND "${GIT_FOUND}" PARENT_SCOPE) + execute_process(COMMAND "${GIT_EXECUTABLE}" describe --abbrev=0 + WORKING_DIRECTORY "${CMAKE_SOURCE_DIR}" + RESULT_VARIABLE _git_status + OUTPUT_VARIABLE _git_output + OUTPUT_STRIP_TRAILING_WHITESPACE) + if((_git_status STREQUAL "0") AND (_git_output MATCHES ${_VERSION_REGEX})) + set(_package_version "${_git_output}") + endif() + execute_process(COMMAND "${GIT_EXECUTABLE}" describe --always + WORKING_DIRECTORY "${CMAKE_SOURCE_DIR}" + RESULT_VARIABLE _git_status + OUTPUT_VARIABLE _full_git_describe + OUTPUT_STRIP_TRAILING_WHITESPACE) + if(NOT (_git_status STREQUAL "0")) + set(_full_git_describe NOTFOUND) + endif() + else() + message( WARNING "Could not find git executable!") + endif() + endif() + + if(NOT (_package_version MATCHES ${_VERSION_REGEX})) + message( WARNING "Could not extract version from git, falling back on ${_VERSION_FILE}.") + file(STRINGS ".VERSION" _package_version + REGEX ${_VERSION_REGEX} + ) + endif() + + if(NOT _full_git_describe) + set(_full_git_describe ${_package_version}) + endif() + + # Strip leading "v" character from package version tags so that + # the version string can be passed to the CMake `project` command + string(REPLACE "v" "" _package_version "${_package_version}") + string(REPLACE "V" "" _package_version "${_package_version}") + + if(set_version_VERSION_VARIABLE) + set(${set_version_VERSION_VARIABLE} ${_package_version} PARENT_SCOPE) + else() + set(PROJECT_VERSION ${_package_version} PARENT_SCOPE) + endif() + if(set_version_GIT_DESCRIBE_VAR) + set(${set_version_GIT_DESCRIBE_VAR} ${_full_git_describe} PARENT_SCOPE) + else() + set(FULL_GIT_DESCRIBE ${_full_git_describe} PARENT_SCOPE) + endif() + +endfunction() diff --git a/cmake/uninstall.cmake.in b/cmake/uninstall.cmake.in new file mode 100644 index 000000000..dd3959305 --- /dev/null +++ b/cmake/uninstall.cmake.in @@ -0,0 +1,23 @@ +# Adapted from http://www.cmake.org/Wiki/CMake_FAQ#Can_I_do_.22make_uninstall.22_with_CMake.3F May 1, 2014 + +if(NOT EXISTS "@CMAKE_BINARY_DIR@/install_manifest.txt") + message(FATAL_ERROR "Cannot find install manifest: @CMAKE_BINARY_DIR@/install_manifest.txt") +endif(NOT EXISTS "@CMAKE_BINARY_DIR@/install_manifest.txt") + +file(READ "@CMAKE_BINARY_DIR@/install_manifest.txt" files) +string(REGEX REPLACE "\n" ";" files "${files}") +foreach(file ${files}) + message(STATUS "Uninstalling $ENV{DESTDIR}${file}") + if(IS_SYMLINK "$ENV{DESTDIR}${file}" OR EXISTS "$ENV{DESTDIR}${file}") + exec_program( + "@CMAKE_COMMAND@" ARGS "-E remove \"$ENV{DESTDIR}${file}\"" + OUTPUT_VARIABLE rm_out + RETURN_VALUE rm_retval + ) + if(NOT "${rm_retval}" STREQUAL 0) + message(FATAL_ERROR "Problem when removing $ENV{DESTDIR}${file}") + endif(NOT "${rm_retval}" STREQUAL 0) + else(IS_SYMLINK "$ENV{DESTDIR}${file}" OR EXISTS "$ENV{DESTDIR}${file}") + message(STATUS "File $ENV{DESTDIR}${file} does not exist.") + endif(IS_SYMLINK "$ENV{DESTDIR}${file}" OR EXISTS "$ENV{DESTDIR}${file}") +endforeach(file) diff --git a/krylov/CMakeLists.txt b/krylov/CMakeLists.txt new file mode 100644 index 000000000..9f2fe2d5d --- /dev/null +++ b/krylov/CMakeLists.txt @@ -0,0 +1,48 @@ +set(PSB_krylov_source_files + psb_base_krylov_conv_mod.f90 + psb_cbicg.f90 + psb_ccg.F90 + psb_ccgs.f90 + psb_ccgstab.f90 + psb_ccgstabl.f90 + psb_cfcg.F90 + psb_cgcr.f90 + psb_c_krylov_conv_mod.f90 + psb_ckrylov.f90 + psb_crgmres.f90 + psb_dbicg.f90 + psb_dcg.F90 + psb_dcgs.f90 + psb_dcgstab.f90 + psb_dcgstabl.f90 + psb_dfcg.F90 + psb_dgcr.f90 + psb_d_krylov_conv_mod.f90 + psb_dkrylov.f90 + psb_drgmres.f90 + psb_krylov_conv_mod.f90 + psb_krylov_mod.f90 + psb_sbicg.f90 + psb_scg.F90 + psb_scgs.f90 + psb_scgstab.f90 + psb_scgstabl.f90 + psb_sfcg.F90 + psb_sgcr.f90 + psb_s_krylov_conv_mod.f90 + psb_skrylov.f90 + psb_srgmres.f90 + psb_zbicg.f90 + psb_zcg.F90 + psb_zcgs.f90 + psb_zcgstab.f90 + psb_zcgstabl.f90 + psb_zfcg.F90 + psb_zgcr.f90 + psb_z_krylov_conv_mod.f90 + psb_zkrylov.f90 + psb_zrgmres.f90 + ) +foreach(file IN LISTS PSB_krylov_source_files) + list(APPEND krylov_source_files ${CMAKE_CURRENT_LIST_DIR}/${file}) +endforeach() diff --git a/krylov/psb_cgcr.f90 b/krylov/psb_cgcr.f90 index 76f581605..2d537bde7 100644 --- a/krylov/psb_cgcr.f90 +++ b/krylov/psb_cgcr.f90 @@ -1,371 +1,369 @@ -! -! Parallel Sparse BLAS version 3.5 -! (C) Copyright 2006, 2010, 2015, 2017 -! Salvatore Filippone Cranfield University -! Alfredo Buttari CNRS-IRIT, Toulouse -! -! 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 PSBLAS group or the names of its contributors may -! not be used to endorse or promote products derived from this -! software without specific 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 PSBLAS 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: psb_cgcr.f90 -!! -!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) -!! -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! C C -! C References: C -! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C -! C Level 3 basic linear algebra subprograms for sparse C -! C matrices: a user level interface C -! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C -! C C -! C C -! C [2] S. Filippone, M. Colajanni C -! C PSBLAS: A library for parallel linear algebra C -! C computation on sparse matrices C -! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C -! C C -! C [3] M. Arioli, I. Duff, M. Ruiz C -! C Stopping criteria for iterative solvers C -! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C -! C C -! C C -! C [4] R. Barrett et al C -! C Templates for the solution of linear systems C -! C SIAM, 1993 -! C C -! C [4] Notay, Yvan C -! C Aggregation-based algebraic multigrid method C -! C SIAM Journal on Scientific Computing 34, C -! C pp. A2288-A2316, 2012 C -! C C -! C C -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! File: psb_cgcr.f90 -! -! Subroutine: psb_cgcr -! This subroutine implements the GCR method. -! -! -! Arguments: -! -! a - type(psb_cspmat_type) Input: sparse matrix containing A. -! prec - class(psb_cprec_type) Input: preconditioner -! b(:) - real Input: vector containing the -! right hand side B -! x(:) - real Input/Output: vector containing the -! initial guess and final solution X. -! eps - real Input: Stopping tolerance; the iteration is -! stopped when the error estimate |err| <= eps -! desc_a - type(psb_desc_type). Input: The communication descriptor. -! info - integer. Output: Return code -! -! itmax - integer(optional) Input: maximum number of iterations to be -! performed. -! iter - integer(optional) Output: how many iterations have been -! performed. -! performed. -! err - real (optional) Output: error estimate on exit. If the -! denominator of the estimate is exactly -! 0, it is changed into 1. -! itrace - integer(optional) Input: print an informational message -! with the error estimate every itrace -! iterations -! istop - integer(optional) Input: stopping criterion, or how -! to estimate the error. -! 1: err = |r|/(|a||x|+|b|); here the iteration is -! stopped when |r| <= eps * (|a||x|+|b|) -! 2: err = |r|/|b|; here the iteration is -! stopped when |r| <= eps * |b| -! where r is the (preconditioned, recursive -! estimate of) residual. -! -! irst - integer(optional) Input: restart parameter -! - -subroutine psb_cgcr_vect(a,prec,b,x,eps,desc_a,info,& - & itmax,iter,err,itrace, irst, istop) - use psb_base_mod - use psb_prec_mod - use psb_c_krylov_conv_mod - use psb_krylov_mod - implicit none - - - type(psb_cspmat_type), intent(in) :: a - Type(psb_desc_type), Intent(in) :: desc_a - class(psb_cprec_type), intent(inout) :: prec - type(psb_c_vect_type), Intent(inout) :: b - type(psb_c_vect_type), Intent(inout) :: x - real(psb_spk_), Intent(in) :: eps - integer(psb_ipk_), intent(out) :: info - integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop - integer(psb_ipk_), Optional, Intent(out) :: iter - real(psb_spk_), Optional, Intent(out) :: err - ! = local data - complex(psb_spk_), allocatable :: alpha(:), h(:,:) - type(psb_c_vect_type), allocatable :: z(:), c(:), c_scale(:) - type(psb_c_vect_type) :: r - - real(psb_dpk_) :: r_norm, b_norm, a_norm, derr - integer(psb_ipk_) :: n_col, mglob, naux, err_act - integer(psb_ipk_) :: debug_level, debug_unit - integer(psb_ipk_) :: np, me, ictxt - integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst - complex(psb_spk_) :: hjj - complex(psb_spk_), allocatable, target :: aux(:) - character(len=20) :: name - type(psb_itconv_type) :: stopdat - character(len=*), parameter :: methdname='GCR' - integer(psb_ipk_) ::int_err(5) - info = psb_success_ - name = 'psb_cgcr' - call psb_erractionsave(err_act) - debug_unit = psb_get_debug_unit() - debug_level = psb_get_debug_level() - - ictxt = desc_a%get_context() - - call psb_info(ictxt, me, np) - if (.not.allocated(b%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - if (.not.allocated(x%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - - mglob = desc_a%get_global_rows() - n_col = desc_a%get_local_cols() - - if (present(istop)) then - istop_ = istop - else - istop_ = 2 - endif - - ! - ! ISTOP_ = 1: Normwise backward error, infinity norm - ! ISTOP_ = 2: ||r||/||b||, 2-norm - ! - - if ((istop_ < 1 ).or.(istop_ > 2 ) ) then - info=psb_err_invalid_istop_ - int_err(1)=istop_ - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - - call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) - if (info == psb_success_)& - & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='psb_chkvect on X/B') - goto 9999 - end if - - if (present(itmax)) then - itmax_ = itmax - else - itmax_ = 1000 - endif - - if (present(itrace)) then - itrace_ = itrace - else - itrace_ = 0 - end if - - if (present(irst)) then - nl = irst - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' present: irst: ',irst,nl - else - nl = 10 - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' not present: irst: ',irst,nl - endif - - if (nl <=0 ) then - info=psb_err_invalid_istop_ - int_err(1)=nl - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - naux=4*n_col - allocate(aux(naux),h(nl+1,nl+1),& - &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) - - h = czero - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) - - do i =1,nl+1 - call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) - end do - - itx = 0 - - nrst = -1 - call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) - restart: do - if (itx>= itmax_) exit restart - h = czero - - it = 0 - ! compute r0 = b-ax0 - - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - - call psb_geaxpby(cone, b, czero, r, desc_a, info) - call psb_spmm(-cone,a,x,cone,r,desc_a,info,work=aux) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - nrst = nrst + 1 - - iteration: do - - itx = itx + 1 - it = it + 1 - j = it - !Apply preconditioner - call prec%apply(r,z(j),desc_a,info,work=aux) - - call psb_spmm(cone,a,z(j),czero,c(1),desc_a,info,work=aux) - do i =1, j - 1 - - h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) - - call psb_geaxpby(cone, c(i), czero, c(i+1), desc_a, info) - call psb_geaxpby(-h(i,j), c_scale(i), cone, c(i+1), desc_a, info) - end do - - h(j,j) = psb_norm2(c(j), desc_a, info) - hjj = cone/h(j,j) - call psb_geaxpby(hjj, c(j), czero, c_scale(j), desc_a, info) - - alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) - - !Update residual - call psb_geaxpby(cone, r, czero, r, desc_a, info) - call psb_geaxpby(-alpha(j), c_scale(j), cone, r, desc_a, info) - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - if (j >= irst) exit iteration - - - end do iteration - - m = j - - !Compute solution - - call ctrsm('l','u','n','n',m,1,cone,h,size(h,1),alpha,size(alpha,1)) - - if (nrst == 0 ) then - call x%set(czero) - endif - do i=1,m - call psb_geaxpby(alpha(i), z(i), cone, x, desc_a, info) - enddo - - - - - end do restart - m = j - !Compute solution - call ctrsm('l','u','n','n',m,1,cone,h,size(h,1),alpha,size(alpha,1)) - call x%set(czero) - do i=1,m - call psb_geaxpby(alpha(i), z(i), cone, x, desc_a, info) - enddo - - iter = j - - call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) - if (present(err)) err = derr - - if (info == psb_success_) call psb_gefree(r,desc_a,info) - - do j = 1,m - if (info == psb_success_) call psb_gefree(z(j),desc_a,info) - if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) - enddo - - do i =1,nl+1 - if (info == psb_success_) call psb_gefree(c(i),desc_a,info) - end do - - if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - if (err_act.eq.psb_act_abort_) then - call psb_error() - return - end if - - return -end subroutine psb_cgcr_vect - - +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006, 2010, 2015, 2017 +! Salvatore Filippone Cranfield University +! Alfredo Buttari CNRS-IRIT, Toulouse +! +! 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 PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific 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 PSBLAS 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: psb_cgcr.f90 +!! +!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) +!! +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! C C +! C References: C +! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C +! C Level 3 basic linear algebra subprograms for sparse C +! C matrices: a user level interface C +! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C +! C C +! C C +! C [2] S. Filippone, M. Colajanni C +! C PSBLAS: A library for parallel linear algebra C +! C computation on sparse matrices C +! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C +! C C +! C [3] M. Arioli, I. Duff, M. Ruiz C +! C Stopping criteria for iterative solvers C +! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C +! C C +! C C +! C [4] R. Barrett et al C +! C Templates for the solution of linear systems C +! C SIAM, 1993 +! C C +! C [4] Notay, Yvan C +! C Aggregation-based algebraic multigrid method C +! C SIAM Journal on Scientific Computing 34, C +! C pp. A2288-A2316, 2012 C +! C C +! C C +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! File: psb_cgcr.f90 +! +! Subroutine: psb_cgcr +! This subroutine implements the GCR method. +! +! +! Arguments: +! +! a - type(psb_cspmat_type) Input: sparse matrix containing A. +! prec - class(psb_cprec_type) Input: preconditioner +! b(:) - real Input: vector containing the +! right hand side B +! x(:) - real Input/Output: vector containing the +! initial guess and final solution X. +! eps - real Input: Stopping tolerance; the iteration is +! stopped when the error estimate |err| <= eps +! desc_a - type(psb_desc_type). Input: The communication descriptor. +! info - integer. Output: Return code +! +! itmax - integer(optional) Input: maximum number of iterations to be +! performed. +! iter - integer(optional) Output: how many iterations have been +! performed. +! performed. +! err - real (optional) Output: error estimate on exit. If the +! denominator of the estimate is exactly +! 0, it is changed into 1. +! itrace - integer(optional) Input: print an informational message +! with the error estimate every itrace +! iterations +! istop - integer(optional) Input: stopping criterion, or how +! to estimate the error. +! 1: err = |r|/(|a||x|+|b|); here the iteration is +! stopped when |r| <= eps * (|a||x|+|b|) +! 2: err = |r|/|b|; here the iteration is +! stopped when |r| <= eps * |b| +! where r is the (preconditioned, recursive +! estimate of) residual. +! +! irst - integer(optional) Input: restart parameter +! + +subroutine psb_cgcr_vect(a,prec,b,x,eps,desc_a,info,& + & itmax,iter,err,itrace, irst, istop) + use psb_base_mod + use psb_prec_mod + use psb_c_krylov_conv_mod + use psb_krylov_mod + implicit none + + + type(psb_cspmat_type), intent(in) :: a + Type(psb_desc_type), Intent(in) :: desc_a + class(psb_cprec_type), intent(inout) :: prec + type(psb_c_vect_type), Intent(inout) :: b + type(psb_c_vect_type), Intent(inout) :: x + real(psb_spk_), Intent(in) :: eps + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop + integer(psb_ipk_), Optional, Intent(out) :: iter + real(psb_spk_), Optional, Intent(out) :: err + ! = local data + complex(psb_spk_), allocatable :: alpha(:), h(:,:) + type(psb_c_vect_type), allocatable :: z(:), c(:), c_scale(:) + type(psb_c_vect_type) :: r + + real(psb_dpk_) :: r_norm, b_norm, a_norm, derr + integer(psb_ipk_) :: n_col, mglob, naux, err_act + integer(psb_ipk_) :: debug_level, debug_unit + integer(psb_ipk_) :: np, me, ictxt + integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst + complex(psb_spk_) :: hjj + complex(psb_spk_), allocatable, target :: aux(:) + character(len=20) :: name + type(psb_itconv_type) :: stopdat + character(len=*), parameter :: methdname='GCR' + integer(psb_ipk_) ::int_err(5) + info = psb_success_ + name = 'psb_cgcr' + call psb_erractionsave(err_act) + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + + ictxt = desc_a%get_context() + + call psb_info(ictxt, me, np) + if (.not.allocated(b%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + if (.not.allocated(x%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + + mglob = desc_a%get_global_rows() + n_col = desc_a%get_local_cols() + + if (present(istop)) then + istop_ = istop + else + istop_ = 2 + endif + + ! + ! ISTOP_ = 1: Normwise backward error, infinity norm + ! ISTOP_ = 2: ||r||/||b||, 2-norm + ! + + if ((istop_ < 1 ).or.(istop_ > 2 ) ) then + info=psb_err_invalid_istop_ + int_err(1)=istop_ + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + + call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) + if (info == psb_success_)& + & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) + if(info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='psb_chkvect on X/B') + goto 9999 + end if + + if (present(itmax)) then + itmax_ = itmax + else + itmax_ = 1000 + endif + + if (present(itrace)) then + itrace_ = itrace + else + itrace_ = 0 + end if + + if (present(irst)) then + nl = irst + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' present: irst: ',irst,nl + else + nl = 10 + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' not present: irst: ',irst,nl + endif + + if (nl <=0 ) then + info=psb_err_invalid_istop_ + int_err(1)=nl + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + naux=4*n_col + allocate(aux(naux),h(nl+1,nl+1),& + &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) + + h = czero + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) + + do i =1,nl+1 + call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) + end do + + itx = 0 + + nrst = -1 + call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) + restart: do + if (itx>= itmax_) exit restart + h = czero + + it = 0 + ! compute r0 = b-ax0 + + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + + call psb_geaxpby(cone, b, czero, r, desc_a, info) + call psb_spmm(-cone,a,x,cone,r,desc_a,info,work=aux) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + nrst = nrst + 1 + + iteration: do + + itx = itx + 1 + it = it + 1 + j = it + !Apply preconditioner + call prec%apply(r,z(j),desc_a,info,work=aux) + + call psb_spmm(cone,a,z(j),czero,c(1),desc_a,info,work=aux) + do i =1, j - 1 + + h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) + + call psb_geaxpby(cone, c(i), czero, c(i+1), desc_a, info) + call psb_geaxpby(-h(i,j), c_scale(i), cone, c(i+1), desc_a, info) + end do + + h(j,j) = psb_norm2(c(j), desc_a, info) + hjj = cone/h(j,j) + call psb_geaxpby(hjj, c(j), czero, c_scale(j), desc_a, info) + + alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) + + !Update residual + call psb_geaxpby(cone, r, czero, r, desc_a, info) + call psb_geaxpby(-alpha(j), c_scale(j), cone, r, desc_a, info) + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + if (j >= irst) exit iteration + + + end do iteration + + m = j + + !Compute solution + + call ctrsm('l','u','n','n',m,1,cone,h,size(h,1),alpha,size(alpha,1)) + + if (nrst == 0 ) then + call x%set(czero) + endif + do i=1,m + call psb_geaxpby(alpha(i), z(i), cone, x, desc_a, info) + enddo + + + + + end do restart + m = j + !Compute solution + call ctrsm('l','u','n','n',m,1,cone,h,size(h,1),alpha,size(alpha,1)) + call x%set(czero) + do i=1,m + call psb_geaxpby(alpha(i), z(i), cone, x, desc_a, info) + enddo + + iter = j + + call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) + if (present(err)) err = derr + + if (info == psb_success_) call psb_gefree(r,desc_a,info) + + do j = 1,m + if (info == psb_success_) call psb_gefree(z(j),desc_a,info) + if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) + enddo + + do i =1,nl+1 + if (info == psb_success_) call psb_gefree(c(i),desc_a,info) + end do + + if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_erractionrestore(err_act) + return + +9999 continue + call psb_erractionrestore(err_act) + if (err_act.eq.psb_act_abort_) then + call psb_error() + return + end if + + return +end subroutine psb_cgcr_vect diff --git a/krylov/psb_dgcr.f90 b/krylov/psb_dgcr.f90 index 7d90720a1..365684a31 100644 --- a/krylov/psb_dgcr.f90 +++ b/krylov/psb_dgcr.f90 @@ -1,371 +1,369 @@ -! -! Parallel Sparse BLAS version 3.5 -! (C) Copyright 2006, 2010, 2015, 2017 -! Salvatore Filippone Cranfield University -! Alfredo Buttari CNRS-IRIT, Toulouse -! -! 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 PSBLAS group or the names of its contributors may -! not be used to endorse or promote products derived from this -! software without specific 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 PSBLAS 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: psb_dgcr.f90 -!! -!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) -!! -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! C C -! C References: C -! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C -! C Level 3 basic linear algebra subprograms for sparse C -! C matrices: a user level interface C -! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C -! C C -! C C -! C [2] S. Filippone, M. Colajanni C -! C PSBLAS: A library for parallel linear algebra C -! C computation on sparse matrices C -! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C -! C C -! C [3] M. Arioli, I. Duff, M. Ruiz C -! C Stopping criteria for iterative solvers C -! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C -! C C -! C C -! C [4] R. Barrett et al C -! C Templates for the solution of linear systems C -! C SIAM, 1993 -! C C -! C [4] Notay, Yvan C -! C Aggregation-based algebraic multigrid method C -! C SIAM Journal on Scientific Computing 34, C -! C pp. A2288-A2316, 2012 C -! C C -! C C -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! File: psb_dgcr.f90 -! -! Subroutine: psb_dgcr -! This subroutine implements the GCR method. -! -! -! Arguments: -! -! a - type(psb_dspmat_type) Input: sparse matrix containing A. -! prec - class(psb_dprec_type) Input: preconditioner -! b(:) - real Input: vector containing the -! right hand side B -! x(:) - real Input/Output: vector containing the -! initial guess and final solution X. -! eps - real Input: Stopping tolerance; the iteration is -! stopped when the error estimate |err| <= eps -! desc_a - type(psb_desc_type). Input: The communication descriptor. -! info - integer. Output: Return code -! -! itmax - integer(optional) Input: maximum number of iterations to be -! performed. -! iter - integer(optional) Output: how many iterations have been -! performed. -! performed. -! err - real (optional) Output: error estimate on exit. If the -! denominator of the estimate is exactly -! 0, it is changed into 1. -! itrace - integer(optional) Input: print an informational message -! with the error estimate every itrace -! iterations -! istop - integer(optional) Input: stopping criterion, or how -! to estimate the error. -! 1: err = |r|/(|a||x|+|b|); here the iteration is -! stopped when |r| <= eps * (|a||x|+|b|) -! 2: err = |r|/|b|; here the iteration is -! stopped when |r| <= eps * |b| -! where r is the (preconditioned, recursive -! estimate of) residual. -! -! irst - integer(optional) Input: restart parameter -! - -subroutine psb_dgcr_vect(a,prec,b,x,eps,desc_a,info,& - & itmax,iter,err,itrace, irst, istop) - use psb_base_mod - use psb_prec_mod - use psb_d_krylov_conv_mod - use psb_krylov_mod - implicit none - - - type(psb_dspmat_type), intent(in) :: a - Type(psb_desc_type), Intent(in) :: desc_a - class(psb_dprec_type), intent(inout) :: prec - type(psb_d_vect_type), Intent(inout) :: b - type(psb_d_vect_type), Intent(inout) :: x - real(psb_dpk_), Intent(in) :: eps - integer(psb_ipk_), intent(out) :: info - integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop - integer(psb_ipk_), Optional, Intent(out) :: iter - real(psb_dpk_), Optional, Intent(out) :: err - ! = local data - real(psb_dpk_), allocatable :: alpha(:), h(:,:) - type(psb_d_vect_type), allocatable :: z(:), c(:), c_scale(:) - type(psb_d_vect_type) :: r - - real(psb_dpk_) :: r_norm, b_norm, a_norm, derr - integer(psb_ipk_) :: n_col, mglob, naux, err_act - integer(psb_ipk_) :: debug_level, debug_unit - integer(psb_ipk_) :: np, me, ictxt - integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst - real(psb_dpk_) :: hjj - real(psb_dpk_), allocatable, target :: aux(:) - character(len=20) :: name - type(psb_itconv_type) :: stopdat - character(len=*), parameter :: methdname='GCR' - integer(psb_ipk_) ::int_err(5) - info = psb_success_ - name = 'psb_dgcr' - call psb_erractionsave(err_act) - debug_unit = psb_get_debug_unit() - debug_level = psb_get_debug_level() - - ictxt = desc_a%get_context() - - call psb_info(ictxt, me, np) - if (.not.allocated(b%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - if (.not.allocated(x%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - - mglob = desc_a%get_global_rows() - n_col = desc_a%get_local_cols() - - if (present(istop)) then - istop_ = istop - else - istop_ = 2 - endif - - ! - ! ISTOP_ = 1: Normwise backward error, infinity norm - ! ISTOP_ = 2: ||r||/||b||, 2-norm - ! - - if ((istop_ < 1 ).or.(istop_ > 2 ) ) then - info=psb_err_invalid_istop_ - int_err(1)=istop_ - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - - call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) - if (info == psb_success_)& - & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='psb_chkvect on X/B') - goto 9999 - end if - - if (present(itmax)) then - itmax_ = itmax - else - itmax_ = 1000 - endif - - if (present(itrace)) then - itrace_ = itrace - else - itrace_ = 0 - end if - - if (present(irst)) then - nl = irst - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' present: irst: ',irst,nl - else - nl = 10 - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' not present: irst: ',irst,nl - endif - - if (nl <=0 ) then - info=psb_err_invalid_istop_ - int_err(1)=nl - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - naux=4*n_col - allocate(aux(naux),h(nl+1,nl+1),& - &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) - - h = dzero - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) - - do i =1,nl+1 - call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) - end do - - itx = 0 - - nrst = -1 - call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) - restart: do - if (itx>= itmax_) exit restart - h = dzero - - it = 0 - ! compute r0 = b-ax0 - - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - - call psb_geaxpby(done, b, dzero, r, desc_a, info) - call psb_spmm(-done,a,x,done,r,desc_a,info,work=aux) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - nrst = nrst + 1 - - iteration: do - - itx = itx + 1 - it = it + 1 - j = it - !Apply preconditioner - call prec%apply(r,z(j),desc_a,info,work=aux) - - call psb_spmm(done,a,z(j),dzero,c(1),desc_a,info,work=aux) - do i =1, j - 1 - - h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) - - call psb_geaxpby(done, c(i), dzero, c(i+1), desc_a, info) - call psb_geaxpby(-h(i,j), c_scale(i), done, c(i+1), desc_a, info) - end do - - h(j,j) = psb_norm2(c(j), desc_a, info) - hjj = done/h(j,j) - call psb_geaxpby(hjj, c(j), dzero, c_scale(j), desc_a, info) - - alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) - - !Update residual - call psb_geaxpby(done, r, dzero, r, desc_a, info) - call psb_geaxpby(-alpha(j), c_scale(j), done, r, desc_a, info) - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - if (j >= irst) exit iteration - - - end do iteration - - m = j - - !Compute solution - - call dtrsm('l','u','n','n',m,1,done,h,size(h,1),alpha,size(alpha,1)) - - if (nrst == 0 ) then - call x%set(dzero) - endif - do i=1,m - call psb_geaxpby(alpha(i), z(i), done, x, desc_a, info) - enddo - - - - - end do restart - m = j - !Compute solution - call dtrsm('l','u','n','n',m,1,done,h,size(h,1),alpha,size(alpha,1)) - call x%set(dzero) - do i=1,m - call psb_geaxpby(alpha(i), z(i), done, x, desc_a, info) - enddo - - iter = j - - call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) - if (present(err)) err = derr - - if (info == psb_success_) call psb_gefree(r,desc_a,info) - - do j = 1,m - if (info == psb_success_) call psb_gefree(z(j),desc_a,info) - if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) - enddo - - do i =1,nl+1 - if (info == psb_success_) call psb_gefree(c(i),desc_a,info) - end do - - if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - if (err_act.eq.psb_act_abort_) then - call psb_error() - return - end if - - return -end subroutine psb_dgcr_vect - - +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006, 2010, 2015, 2017 +! Salvatore Filippone Cranfield University +! Alfredo Buttari CNRS-IRIT, Toulouse +! +! 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 PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific 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 PSBLAS 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: psb_dgcr.f90 +!! +!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) +!! +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! C C +! C References: C +! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C +! C Level 3 basic linear algebra subprograms for sparse C +! C matrices: a user level interface C +! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C +! C C +! C C +! C [2] S. Filippone, M. Colajanni C +! C PSBLAS: A library for parallel linear algebra C +! C computation on sparse matrices C +! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C +! C C +! C [3] M. Arioli, I. Duff, M. Ruiz C +! C Stopping criteria for iterative solvers C +! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C +! C C +! C C +! C [4] R. Barrett et al C +! C Templates for the solution of linear systems C +! C SIAM, 1993 +! C C +! C [4] Notay, Yvan C +! C Aggregation-based algebraic multigrid method C +! C SIAM Journal on Scientific Computing 34, C +! C pp. A2288-A2316, 2012 C +! C C +! C C +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! File: psb_dgcr.f90 +! +! Subroutine: psb_dgcr +! This subroutine implements the GCR method. +! +! +! Arguments: +! +! a - type(psb_dspmat_type) Input: sparse matrix containing A. +! prec - class(psb_dprec_type) Input: preconditioner +! b(:) - real Input: vector containing the +! right hand side B +! x(:) - real Input/Output: vector containing the +! initial guess and final solution X. +! eps - real Input: Stopping tolerance; the iteration is +! stopped when the error estimate |err| <= eps +! desc_a - type(psb_desc_type). Input: The communication descriptor. +! info - integer. Output: Return code +! +! itmax - integer(optional) Input: maximum number of iterations to be +! performed. +! iter - integer(optional) Output: how many iterations have been +! performed. +! performed. +! err - real (optional) Output: error estimate on exit. If the +! denominator of the estimate is exactly +! 0, it is changed into 1. +! itrace - integer(optional) Input: print an informational message +! with the error estimate every itrace +! iterations +! istop - integer(optional) Input: stopping criterion, or how +! to estimate the error. +! 1: err = |r|/(|a||x|+|b|); here the iteration is +! stopped when |r| <= eps * (|a||x|+|b|) +! 2: err = |r|/|b|; here the iteration is +! stopped when |r| <= eps * |b| +! where r is the (preconditioned, recursive +! estimate of) residual. +! +! irst - integer(optional) Input: restart parameter +! + +subroutine psb_dgcr_vect(a,prec,b,x,eps,desc_a,info,& + & itmax,iter,err,itrace, irst, istop) + use psb_base_mod + use psb_prec_mod + use psb_d_krylov_conv_mod + use psb_krylov_mod + implicit none + + + type(psb_dspmat_type), intent(in) :: a + Type(psb_desc_type), Intent(in) :: desc_a + class(psb_dprec_type), intent(inout) :: prec + type(psb_d_vect_type), Intent(inout) :: b + type(psb_d_vect_type), Intent(inout) :: x + real(psb_dpk_), Intent(in) :: eps + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop + integer(psb_ipk_), Optional, Intent(out) :: iter + real(psb_dpk_), Optional, Intent(out) :: err + ! = local data + real(psb_dpk_), allocatable :: alpha(:), h(:,:) + type(psb_d_vect_type), allocatable :: z(:), c(:), c_scale(:) + type(psb_d_vect_type) :: r + + real(psb_dpk_) :: r_norm, b_norm, a_norm, derr + integer(psb_ipk_) :: n_col, mglob, naux, err_act + integer(psb_ipk_) :: debug_level, debug_unit + integer(psb_ipk_) :: np, me, ictxt + integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst + real(psb_dpk_) :: hjj + real(psb_dpk_), allocatable, target :: aux(:) + character(len=20) :: name + type(psb_itconv_type) :: stopdat + character(len=*), parameter :: methdname='GCR' + integer(psb_ipk_) ::int_err(5) + info = psb_success_ + name = 'psb_dgcr' + call psb_erractionsave(err_act) + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + + ictxt = desc_a%get_context() + + call psb_info(ictxt, me, np) + if (.not.allocated(b%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + if (.not.allocated(x%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + + mglob = desc_a%get_global_rows() + n_col = desc_a%get_local_cols() + + if (present(istop)) then + istop_ = istop + else + istop_ = 2 + endif + + ! + ! ISTOP_ = 1: Normwise backward error, infinity norm + ! ISTOP_ = 2: ||r||/||b||, 2-norm + ! + + if ((istop_ < 1 ).or.(istop_ > 2 ) ) then + info=psb_err_invalid_istop_ + int_err(1)=istop_ + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + + call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) + if (info == psb_success_)& + & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) + if(info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='psb_chkvect on X/B') + goto 9999 + end if + + if (present(itmax)) then + itmax_ = itmax + else + itmax_ = 1000 + endif + + if (present(itrace)) then + itrace_ = itrace + else + itrace_ = 0 + end if + + if (present(irst)) then + nl = irst + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' present: irst: ',irst,nl + else + nl = 10 + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' not present: irst: ',irst,nl + endif + + if (nl <=0 ) then + info=psb_err_invalid_istop_ + int_err(1)=nl + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + naux=4*n_col + allocate(aux(naux),h(nl+1,nl+1),& + &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) + + h = dzero + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) + + do i =1,nl+1 + call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) + end do + + itx = 0 + + nrst = -1 + call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) + restart: do + if (itx>= itmax_) exit restart + h = dzero + + it = 0 + ! compute r0 = b-ax0 + + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + + call psb_geaxpby(done, b, dzero, r, desc_a, info) + call psb_spmm(-done,a,x,done,r,desc_a,info,work=aux) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + nrst = nrst + 1 + + iteration: do + + itx = itx + 1 + it = it + 1 + j = it + !Apply preconditioner + call prec%apply(r,z(j),desc_a,info,work=aux) + + call psb_spmm(done,a,z(j),dzero,c(1),desc_a,info,work=aux) + do i =1, j - 1 + + h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) + + call psb_geaxpby(done, c(i), dzero, c(i+1), desc_a, info) + call psb_geaxpby(-h(i,j), c_scale(i), done, c(i+1), desc_a, info) + end do + + h(j,j) = psb_norm2(c(j), desc_a, info) + hjj = done/h(j,j) + call psb_geaxpby(hjj, c(j), dzero, c_scale(j), desc_a, info) + + alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) + + !Update residual + call psb_geaxpby(done, r, dzero, r, desc_a, info) + call psb_geaxpby(-alpha(j), c_scale(j), done, r, desc_a, info) + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + if (j >= irst) exit iteration + + + end do iteration + + m = j + + !Compute solution + + call dtrsm('l','u','n','n',m,1,done,h,size(h,1),alpha,size(alpha,1)) + + if (nrst == 0 ) then + call x%set(dzero) + endif + do i=1,m + call psb_geaxpby(alpha(i), z(i), done, x, desc_a, info) + enddo + + + + + end do restart + m = j + !Compute solution + call dtrsm('l','u','n','n',m,1,done,h,size(h,1),alpha,size(alpha,1)) + call x%set(dzero) + do i=1,m + call psb_geaxpby(alpha(i), z(i), done, x, desc_a, info) + enddo + + iter = j + + call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) + if (present(err)) err = derr + + if (info == psb_success_) call psb_gefree(r,desc_a,info) + + do j = 1,m + if (info == psb_success_) call psb_gefree(z(j),desc_a,info) + if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) + enddo + + do i =1,nl+1 + if (info == psb_success_) call psb_gefree(c(i),desc_a,info) + end do + + if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_erractionrestore(err_act) + return + +9999 continue + call psb_erractionrestore(err_act) + if (err_act.eq.psb_act_abort_) then + call psb_error() + return + end if + + return +end subroutine psb_dgcr_vect diff --git a/krylov/psb_sgcr.f90 b/krylov/psb_sgcr.f90 index 9bc81f709..23a9326ba 100644 --- a/krylov/psb_sgcr.f90 +++ b/krylov/psb_sgcr.f90 @@ -1,371 +1,369 @@ -! -! Parallel Sparse BLAS version 3.5 -! (C) Copyright 2006, 2010, 2015, 2017 -! Salvatore Filippone Cranfield University -! Alfredo Buttari CNRS-IRIT, Toulouse -! -! 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 PSBLAS group or the names of its contributors may -! not be used to endorse or promote products derived from this -! software without specific 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 PSBLAS 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: psb_sgcr.f90 -!! -!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) -!! -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! C C -! C References: C -! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C -! C Level 3 basic linear algebra subprograms for sparse C -! C matrices: a user level interface C -! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C -! C C -! C C -! C [2] S. Filippone, M. Colajanni C -! C PSBLAS: A library for parallel linear algebra C -! C computation on sparse matrices C -! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C -! C C -! C [3] M. Arioli, I. Duff, M. Ruiz C -! C Stopping criteria for iterative solvers C -! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C -! C C -! C C -! C [4] R. Barrett et al C -! C Templates for the solution of linear systems C -! C SIAM, 1993 -! C C -! C [4] Notay, Yvan C -! C Aggregation-based algebraic multigrid method C -! C SIAM Journal on Scientific Computing 34, C -! C pp. A2288-A2316, 2012 C -! C C -! C C -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! File: psb_sgcr.f90 -! -! Subroutine: psb_sgcr -! This subroutine implements the GCR method. -! -! -! Arguments: -! -! a - type(psb_sspmat_type) Input: sparse matrix containing A. -! prec - class(psb_sprec_type) Input: preconditioner -! b(:) - real Input: vector containing the -! right hand side B -! x(:) - real Input/Output: vector containing the -! initial guess and final solution X. -! eps - real Input: Stopping tolerance; the iteration is -! stopped when the error estimate |err| <= eps -! desc_a - type(psb_desc_type). Input: The communication descriptor. -! info - integer. Output: Return code -! -! itmax - integer(optional) Input: maximum number of iterations to be -! performed. -! iter - integer(optional) Output: how many iterations have been -! performed. -! performed. -! err - real (optional) Output: error estimate on exit. If the -! denominator of the estimate is exactly -! 0, it is changed into 1. -! itrace - integer(optional) Input: print an informational message -! with the error estimate every itrace -! iterations -! istop - integer(optional) Input: stopping criterion, or how -! to estimate the error. -! 1: err = |r|/(|a||x|+|b|); here the iteration is -! stopped when |r| <= eps * (|a||x|+|b|) -! 2: err = |r|/|b|; here the iteration is -! stopped when |r| <= eps * |b| -! where r is the (preconditioned, recursive -! estimate of) residual. -! -! irst - integer(optional) Input: restart parameter -! - -subroutine psb_sgcr_vect(a,prec,b,x,eps,desc_a,info,& - & itmax,iter,err,itrace, irst, istop) - use psb_base_mod - use psb_prec_mod - use psb_s_krylov_conv_mod - use psb_krylov_mod - implicit none - - - type(psb_sspmat_type), intent(in) :: a - Type(psb_desc_type), Intent(in) :: desc_a - class(psb_sprec_type), intent(inout) :: prec - type(psb_s_vect_type), Intent(inout) :: b - type(psb_s_vect_type), Intent(inout) :: x - real(psb_spk_), Intent(in) :: eps - integer(psb_ipk_), intent(out) :: info - integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop - integer(psb_ipk_), Optional, Intent(out) :: iter - real(psb_spk_), Optional, Intent(out) :: err - ! = local data - real(psb_spk_), allocatable :: alpha(:), h(:,:) - type(psb_s_vect_type), allocatable :: z(:), c(:), c_scale(:) - type(psb_s_vect_type) :: r - - real(psb_dpk_) :: r_norm, b_norm, a_norm, derr - integer(psb_ipk_) :: n_col, mglob, naux, err_act - integer(psb_ipk_) :: debug_level, debug_unit - integer(psb_ipk_) :: np, me, ictxt - integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst - real(psb_spk_) :: hjj - real(psb_spk_), allocatable, target :: aux(:) - character(len=20) :: name - type(psb_itconv_type) :: stopdat - character(len=*), parameter :: methdname='GCR' - integer(psb_ipk_) ::int_err(5) - info = psb_success_ - name = 'psb_sgcr' - call psb_erractionsave(err_act) - debug_unit = psb_get_debug_unit() - debug_level = psb_get_debug_level() - - ictxt = desc_a%get_context() - - call psb_info(ictxt, me, np) - if (.not.allocated(b%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - if (.not.allocated(x%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - - mglob = desc_a%get_global_rows() - n_col = desc_a%get_local_cols() - - if (present(istop)) then - istop_ = istop - else - istop_ = 2 - endif - - ! - ! ISTOP_ = 1: Normwise backward error, infinity norm - ! ISTOP_ = 2: ||r||/||b||, 2-norm - ! - - if ((istop_ < 1 ).or.(istop_ > 2 ) ) then - info=psb_err_invalid_istop_ - int_err(1)=istop_ - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - - call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) - if (info == psb_success_)& - & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='psb_chkvect on X/B') - goto 9999 - end if - - if (present(itmax)) then - itmax_ = itmax - else - itmax_ = 1000 - endif - - if (present(itrace)) then - itrace_ = itrace - else - itrace_ = 0 - end if - - if (present(irst)) then - nl = irst - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' present: irst: ',irst,nl - else - nl = 10 - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' not present: irst: ',irst,nl - endif - - if (nl <=0 ) then - info=psb_err_invalid_istop_ - int_err(1)=nl - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - naux=4*n_col - allocate(aux(naux),h(nl+1,nl+1),& - &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) - - h = szero - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) - - do i =1,nl+1 - call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) - end do - - itx = 0 - - nrst = -1 - call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) - restart: do - if (itx>= itmax_) exit restart - h = szero - - it = 0 - ! compute r0 = b-ax0 - - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - - call psb_geaxpby(sone, b, szero, r, desc_a, info) - call psb_spmm(-sone,a,x,sone,r,desc_a,info,work=aux) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - nrst = nrst + 1 - - iteration: do - - itx = itx + 1 - it = it + 1 - j = it - !Apply preconditioner - call prec%apply(r,z(j),desc_a,info,work=aux) - - call psb_spmm(sone,a,z(j),szero,c(1),desc_a,info,work=aux) - do i =1, j - 1 - - h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) - - call psb_geaxpby(sone, c(i), szero, c(i+1), desc_a, info) - call psb_geaxpby(-h(i,j), c_scale(i), sone, c(i+1), desc_a, info) - end do - - h(j,j) = psb_norm2(c(j), desc_a, info) - hjj = sone/h(j,j) - call psb_geaxpby(hjj, c(j), szero, c_scale(j), desc_a, info) - - alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) - - !Update residual - call psb_geaxpby(sone, r, szero, r, desc_a, info) - call psb_geaxpby(-alpha(j), c_scale(j), sone, r, desc_a, info) - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - if (j >= irst) exit iteration - - - end do iteration - - m = j - - !Compute solution - - call strsm('l','u','n','n',m,1,sone,h,size(h,1),alpha,size(alpha,1)) - - if (nrst == 0 ) then - call x%set(szero) - endif - do i=1,m - call psb_geaxpby(alpha(i), z(i), sone, x, desc_a, info) - enddo - - - - - end do restart - m = j - !Compute solution - call strsm('l','u','n','n',m,1,sone,h,size(h,1),alpha,size(alpha,1)) - call x%set(szero) - do i=1,m - call psb_geaxpby(alpha(i), z(i), sone, x, desc_a, info) - enddo - - iter = j - - call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) - if (present(err)) err = derr - - if (info == psb_success_) call psb_gefree(r,desc_a,info) - - do j = 1,m - if (info == psb_success_) call psb_gefree(z(j),desc_a,info) - if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) - enddo - - do i =1,nl+1 - if (info == psb_success_) call psb_gefree(c(i),desc_a,info) - end do - - if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - if (err_act.eq.psb_act_abort_) then - call psb_error() - return - end if - - return -end subroutine psb_sgcr_vect - - +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006, 2010, 2015, 2017 +! Salvatore Filippone Cranfield University +! Alfredo Buttari CNRS-IRIT, Toulouse +! +! 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 PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific 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 PSBLAS 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: psb_sgcr.f90 +!! +!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) +!! +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! C C +! C References: C +! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C +! C Level 3 basic linear algebra subprograms for sparse C +! C matrices: a user level interface C +! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C +! C C +! C C +! C [2] S. Filippone, M. Colajanni C +! C PSBLAS: A library for parallel linear algebra C +! C computation on sparse matrices C +! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C +! C C +! C [3] M. Arioli, I. Duff, M. Ruiz C +! C Stopping criteria for iterative solvers C +! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C +! C C +! C C +! C [4] R. Barrett et al C +! C Templates for the solution of linear systems C +! C SIAM, 1993 +! C C +! C [4] Notay, Yvan C +! C Aggregation-based algebraic multigrid method C +! C SIAM Journal on Scientific Computing 34, C +! C pp. A2288-A2316, 2012 C +! C C +! C C +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! File: psb_sgcr.f90 +! +! Subroutine: psb_sgcr +! This subroutine implements the GCR method. +! +! +! Arguments: +! +! a - type(psb_sspmat_type) Input: sparse matrix containing A. +! prec - class(psb_sprec_type) Input: preconditioner +! b(:) - real Input: vector containing the +! right hand side B +! x(:) - real Input/Output: vector containing the +! initial guess and final solution X. +! eps - real Input: Stopping tolerance; the iteration is +! stopped when the error estimate |err| <= eps +! desc_a - type(psb_desc_type). Input: The communication descriptor. +! info - integer. Output: Return code +! +! itmax - integer(optional) Input: maximum number of iterations to be +! performed. +! iter - integer(optional) Output: how many iterations have been +! performed. +! performed. +! err - real (optional) Output: error estimate on exit. If the +! denominator of the estimate is exactly +! 0, it is changed into 1. +! itrace - integer(optional) Input: print an informational message +! with the error estimate every itrace +! iterations +! istop - integer(optional) Input: stopping criterion, or how +! to estimate the error. +! 1: err = |r|/(|a||x|+|b|); here the iteration is +! stopped when |r| <= eps * (|a||x|+|b|) +! 2: err = |r|/|b|; here the iteration is +! stopped when |r| <= eps * |b| +! where r is the (preconditioned, recursive +! estimate of) residual. +! +! irst - integer(optional) Input: restart parameter +! + +subroutine psb_sgcr_vect(a,prec,b,x,eps,desc_a,info,& + & itmax,iter,err,itrace, irst, istop) + use psb_base_mod + use psb_prec_mod + use psb_s_krylov_conv_mod + use psb_krylov_mod + implicit none + + + type(psb_sspmat_type), intent(in) :: a + Type(psb_desc_type), Intent(in) :: desc_a + class(psb_sprec_type), intent(inout) :: prec + type(psb_s_vect_type), Intent(inout) :: b + type(psb_s_vect_type), Intent(inout) :: x + real(psb_spk_), Intent(in) :: eps + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop + integer(psb_ipk_), Optional, Intent(out) :: iter + real(psb_spk_), Optional, Intent(out) :: err + ! = local data + real(psb_spk_), allocatable :: alpha(:), h(:,:) + type(psb_s_vect_type), allocatable :: z(:), c(:), c_scale(:) + type(psb_s_vect_type) :: r + + real(psb_dpk_) :: r_norm, b_norm, a_norm, derr + integer(psb_ipk_) :: n_col, mglob, naux, err_act + integer(psb_ipk_) :: debug_level, debug_unit + integer(psb_ipk_) :: np, me, ictxt + integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst + real(psb_spk_) :: hjj + real(psb_spk_), allocatable, target :: aux(:) + character(len=20) :: name + type(psb_itconv_type) :: stopdat + character(len=*), parameter :: methdname='GCR' + integer(psb_ipk_) ::int_err(5) + info = psb_success_ + name = 'psb_sgcr' + call psb_erractionsave(err_act) + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + + ictxt = desc_a%get_context() + + call psb_info(ictxt, me, np) + if (.not.allocated(b%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + if (.not.allocated(x%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + + mglob = desc_a%get_global_rows() + n_col = desc_a%get_local_cols() + + if (present(istop)) then + istop_ = istop + else + istop_ = 2 + endif + + ! + ! ISTOP_ = 1: Normwise backward error, infinity norm + ! ISTOP_ = 2: ||r||/||b||, 2-norm + ! + + if ((istop_ < 1 ).or.(istop_ > 2 ) ) then + info=psb_err_invalid_istop_ + int_err(1)=istop_ + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + + call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) + if (info == psb_success_)& + & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) + if(info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='psb_chkvect on X/B') + goto 9999 + end if + + if (present(itmax)) then + itmax_ = itmax + else + itmax_ = 1000 + endif + + if (present(itrace)) then + itrace_ = itrace + else + itrace_ = 0 + end if + + if (present(irst)) then + nl = irst + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' present: irst: ',irst,nl + else + nl = 10 + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' not present: irst: ',irst,nl + endif + + if (nl <=0 ) then + info=psb_err_invalid_istop_ + int_err(1)=nl + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + naux=4*n_col + allocate(aux(naux),h(nl+1,nl+1),& + &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) + + h = szero + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) + + do i =1,nl+1 + call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) + end do + + itx = 0 + + nrst = -1 + call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) + restart: do + if (itx>= itmax_) exit restart + h = szero + + it = 0 + ! compute r0 = b-ax0 + + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + + call psb_geaxpby(sone, b, szero, r, desc_a, info) + call psb_spmm(-sone,a,x,sone,r,desc_a,info,work=aux) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + nrst = nrst + 1 + + iteration: do + + itx = itx + 1 + it = it + 1 + j = it + !Apply preconditioner + call prec%apply(r,z(j),desc_a,info,work=aux) + + call psb_spmm(sone,a,z(j),szero,c(1),desc_a,info,work=aux) + do i =1, j - 1 + + h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) + + call psb_geaxpby(sone, c(i), szero, c(i+1), desc_a, info) + call psb_geaxpby(-h(i,j), c_scale(i), sone, c(i+1), desc_a, info) + end do + + h(j,j) = psb_norm2(c(j), desc_a, info) + hjj = sone/h(j,j) + call psb_geaxpby(hjj, c(j), szero, c_scale(j), desc_a, info) + + alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) + + !Update residual + call psb_geaxpby(sone, r, szero, r, desc_a, info) + call psb_geaxpby(-alpha(j), c_scale(j), sone, r, desc_a, info) + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + if (j >= irst) exit iteration + + + end do iteration + + m = j + + !Compute solution + + call strsm('l','u','n','n',m,1,sone,h,size(h,1),alpha,size(alpha,1)) + + if (nrst == 0 ) then + call x%set(szero) + endif + do i=1,m + call psb_geaxpby(alpha(i), z(i), sone, x, desc_a, info) + enddo + + + + + end do restart + m = j + !Compute solution + call strsm('l','u','n','n',m,1,sone,h,size(h,1),alpha,size(alpha,1)) + call x%set(szero) + do i=1,m + call psb_geaxpby(alpha(i), z(i), sone, x, desc_a, info) + enddo + + iter = j + + call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) + if (present(err)) err = derr + + if (info == psb_success_) call psb_gefree(r,desc_a,info) + + do j = 1,m + if (info == psb_success_) call psb_gefree(z(j),desc_a,info) + if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) + enddo + + do i =1,nl+1 + if (info == psb_success_) call psb_gefree(c(i),desc_a,info) + end do + + if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_erractionrestore(err_act) + return + +9999 continue + call psb_erractionrestore(err_act) + if (err_act.eq.psb_act_abort_) then + call psb_error() + return + end if + + return +end subroutine psb_sgcr_vect diff --git a/krylov/psb_zgcr.f90 b/krylov/psb_zgcr.f90 index 6d17db930..42cf68dc8 100644 --- a/krylov/psb_zgcr.f90 +++ b/krylov/psb_zgcr.f90 @@ -1,371 +1,369 @@ -! -! Parallel Sparse BLAS version 3.5 -! (C) Copyright 2006, 2010, 2015, 2017 -! Salvatore Filippone Cranfield University -! Alfredo Buttari CNRS-IRIT, Toulouse -! -! 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 PSBLAS group or the names of its contributors may -! not be used to endorse or promote products derived from this -! software without specific 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 PSBLAS 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: psb_zgcr.f90 -!! -!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) -!! -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! C C -! C References: C -! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C -! C Level 3 basic linear algebra subprograms for sparse C -! C matrices: a user level interface C -! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C -! C C -! C C -! C [2] S. Filippone, M. Colajanni C -! C PSBLAS: A library for parallel linear algebra C -! C computation on sparse matrices C -! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C -! C C -! C [3] M. Arioli, I. Duff, M. Ruiz C -! C Stopping criteria for iterative solvers C -! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C -! C C -! C C -! C [4] R. Barrett et al C -! C Templates for the solution of linear systems C -! C SIAM, 1993 -! C C -! C [4] Notay, Yvan C -! C Aggregation-based algebraic multigrid method C -! C SIAM Journal on Scientific Computing 34, C -! C pp. A2288-A2316, 2012 C -! C C -! C C -! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -! File: psb_zgcr.f90 -! -! Subroutine: psb_zgcr -! This subroutine implements the GCR method. -! -! -! Arguments: -! -! a - type(psb_zspmat_type) Input: sparse matrix containing A. -! prec - class(psb_zprec_type) Input: preconditioner -! b(:) - real Input: vector containing the -! right hand side B -! x(:) - real Input/Output: vector containing the -! initial guess and final solution X. -! eps - real Input: Stopping tolerance; the iteration is -! stopped when the error estimate |err| <= eps -! desc_a - type(psb_desc_type). Input: The communication descriptor. -! info - integer. Output: Return code -! -! itmax - integer(optional) Input: maximum number of iterations to be -! performed. -! iter - integer(optional) Output: how many iterations have been -! performed. -! performed. -! err - real (optional) Output: error estimate on exit. If the -! denominator of the estimate is exactly -! 0, it is changed into 1. -! itrace - integer(optional) Input: print an informational message -! with the error estimate every itrace -! iterations -! istop - integer(optional) Input: stopping criterion, or how -! to estimate the error. -! 1: err = |r|/(|a||x|+|b|); here the iteration is -! stopped when |r| <= eps * (|a||x|+|b|) -! 2: err = |r|/|b|; here the iteration is -! stopped when |r| <= eps * |b| -! where r is the (preconditioned, recursive -! estimate of) residual. -! -! irst - integer(optional) Input: restart parameter -! - -subroutine psb_zgcr_vect(a,prec,b,x,eps,desc_a,info,& - & itmax,iter,err,itrace, irst, istop) - use psb_base_mod - use psb_prec_mod - use psb_z_krylov_conv_mod - use psb_krylov_mod - implicit none - - - type(psb_zspmat_type), intent(in) :: a - Type(psb_desc_type), Intent(in) :: desc_a - class(psb_zprec_type), intent(inout) :: prec - type(psb_z_vect_type), Intent(inout) :: b - type(psb_z_vect_type), Intent(inout) :: x - real(psb_dpk_), Intent(in) :: eps - integer(psb_ipk_), intent(out) :: info - integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop - integer(psb_ipk_), Optional, Intent(out) :: iter - real(psb_dpk_), Optional, Intent(out) :: err - ! = local data - complex(psb_dpk_), allocatable :: alpha(:), h(:,:) - type(psb_z_vect_type), allocatable :: z(:), c(:), c_scale(:) - type(psb_z_vect_type) :: r - - real(psb_dpk_) :: r_norm, b_norm, a_norm, derr - integer(psb_ipk_) :: n_col, mglob, naux, err_act - integer(psb_ipk_) :: debug_level, debug_unit - integer(psb_ipk_) :: np, me, ictxt - integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst - complex(psb_dpk_) :: hjj - complex(psb_dpk_), allocatable, target :: aux(:) - character(len=20) :: name - type(psb_itconv_type) :: stopdat - character(len=*), parameter :: methdname='GCR' - integer(psb_ipk_) ::int_err(5) - info = psb_success_ - name = 'psb_zgcr' - call psb_erractionsave(err_act) - debug_unit = psb_get_debug_unit() - debug_level = psb_get_debug_level() - - ictxt = desc_a%get_context() - - call psb_info(ictxt, me, np) - if (.not.allocated(b%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - if (.not.allocated(x%v)) then - info = psb_err_invalid_vect_state_ - call psb_errpush(info,name) - goto 9999 - endif - - mglob = desc_a%get_global_rows() - n_col = desc_a%get_local_cols() - - if (present(istop)) then - istop_ = istop - else - istop_ = 2 - endif - - ! - ! ISTOP_ = 1: Normwise backward error, infinity norm - ! ISTOP_ = 2: ||r||/||b||, 2-norm - ! - - if ((istop_ < 1 ).or.(istop_ > 2 ) ) then - info=psb_err_invalid_istop_ - int_err(1)=istop_ - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - - call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) - if (info == psb_success_)& - & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) - if(info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='psb_chkvect on X/B') - goto 9999 - end if - - if (present(itmax)) then - itmax_ = itmax - else - itmax_ = 1000 - endif - - if (present(itrace)) then - itrace_ = itrace - else - itrace_ = 0 - end if - - if (present(irst)) then - nl = irst - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' present: irst: ',irst,nl - else - nl = 10 - if (debug_level >= psb_debug_ext_) & - & write(debug_unit,*) me,' ',trim(name),& - & ' not present: irst: ',irst,nl - endif - - if (nl <=0 ) then - info=psb_err_invalid_istop_ - int_err(1)=nl - err=info - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - naux=4*n_col - allocate(aux(naux),h(nl+1,nl+1),& - &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) - - h = zzero - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) - - do i =1,nl+1 - call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) - call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) - end do - - itx = 0 - - nrst = -1 - call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) - restart: do - if (itx>= itmax_) exit restart - h = zzero - - it = 0 - ! compute r0 = b-ax0 - - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - - call psb_geaxpby(zone, b, zzero, r, desc_a, info) - call psb_spmm(-zone,a,x,zone,r,desc_a,info,work=aux) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - nrst = nrst + 1 - - iteration: do - - itx = itx + 1 - it = it + 1 - j = it - !Apply preconditioner - call prec%apply(r,z(j),desc_a,info,work=aux) - - call psb_spmm(zone,a,z(j),zzero,c(1),desc_a,info,work=aux) - do i =1, j - 1 - - h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) - - call psb_geaxpby(zone, c(i), zzero, c(i+1), desc_a, info) - call psb_geaxpby(-h(i,j), c_scale(i), zone, c(i+1), desc_a, info) - end do - - h(j,j) = psb_norm2(c(j), desc_a, info) - hjj = zone/h(j,j) - call psb_geaxpby(hjj, c(j), zzero, c_scale(j), desc_a, info) - - alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) - - !Update residual - call psb_geaxpby(zone, r, zzero, r, desc_a, info) - call psb_geaxpby(-alpha(j), c_scale(j), zone, r, desc_a, info) - - if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart - - if (j >= irst) exit iteration - - - end do iteration - - m = j - - !Compute solution - - call ztrsm('l','u','n','n',m,1,zone,h,size(h,1),alpha,size(alpha,1)) - - if (nrst == 0 ) then - call x%set(zzero) - endif - do i=1,m - call psb_geaxpby(alpha(i), z(i), zone, x, desc_a, info) - enddo - - - - - end do restart - m = j - !Compute solution - call ztrsm('l','u','n','n',m,1,zone,h,size(h,1),alpha,size(alpha,1)) - call x%set(zzero) - do i=1,m - call psb_geaxpby(alpha(i), z(i), zone, x, desc_a, info) - enddo - - iter = j - - call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) - if (present(err)) err = derr - - if (info == psb_success_) call psb_gefree(r,desc_a,info) - - do j = 1,m - if (info == psb_success_) call psb_gefree(z(j),desc_a,info) - if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) - enddo - - do i =1,nl+1 - if (info == psb_success_) call psb_gefree(c(i),desc_a,info) - end do - - if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_non_ - call psb_errpush(info,name) - goto 9999 - end if - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - if (err_act.eq.psb_act_abort_) then - call psb_error() - return - end if - - return -end subroutine psb_zgcr_vect - - +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006, 2010, 2015, 2017 +! Salvatore Filippone Cranfield University +! Alfredo Buttari CNRS-IRIT, Toulouse +! +! 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 PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific 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 PSBLAS 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: psb_zgcr.f90 +!! +!! Contributors: Ambra Abdullahi (UNITOV) and Pasqua D’Ambra (IAC-CNR) +!! +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! C C +! C References: C +! C [1] Duff, I., Marrone, M., Radicati, G., and Vittoli, C. C +! C Level 3 basic linear algebra subprograms for sparse C +! C matrices: a user level interface C +! C ACM Trans. Math. Softw., 23(3), 379-401, 1997. C +! C C +! C C +! C [2] S. Filippone, M. Colajanni C +! C PSBLAS: A library for parallel linear algebra C +! C computation on sparse matrices C +! C ACM Trans. on Math. Softw., 26(4), 527-550, Dec. 2000. C +! C C +! C [3] M. Arioli, I. Duff, M. Ruiz C +! C Stopping criteria for iterative solvers C +! C SIAM J. Matrix Anal. Appl., Vol. 13, pp. 138-144, 1992 C +! C C +! C C +! C [4] R. Barrett et al C +! C Templates for the solution of linear systems C +! C SIAM, 1993 +! C C +! C [4] Notay, Yvan C +! C Aggregation-based algebraic multigrid method C +! C SIAM Journal on Scientific Computing 34, C +! C pp. A2288-A2316, 2012 C +! C C +! C C +! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC +! File: psb_zgcr.f90 +! +! Subroutine: psb_zgcr +! This subroutine implements the GCR method. +! +! +! Arguments: +! +! a - type(psb_zspmat_type) Input: sparse matrix containing A. +! prec - class(psb_zprec_type) Input: preconditioner +! b(:) - real Input: vector containing the +! right hand side B +! x(:) - real Input/Output: vector containing the +! initial guess and final solution X. +! eps - real Input: Stopping tolerance; the iteration is +! stopped when the error estimate |err| <= eps +! desc_a - type(psb_desc_type). Input: The communication descriptor. +! info - integer. Output: Return code +! +! itmax - integer(optional) Input: maximum number of iterations to be +! performed. +! iter - integer(optional) Output: how many iterations have been +! performed. +! performed. +! err - real (optional) Output: error estimate on exit. If the +! denominator of the estimate is exactly +! 0, it is changed into 1. +! itrace - integer(optional) Input: print an informational message +! with the error estimate every itrace +! iterations +! istop - integer(optional) Input: stopping criterion, or how +! to estimate the error. +! 1: err = |r|/(|a||x|+|b|); here the iteration is +! stopped when |r| <= eps * (|a||x|+|b|) +! 2: err = |r|/|b|; here the iteration is +! stopped when |r| <= eps * |b| +! where r is the (preconditioned, recursive +! estimate of) residual. +! +! irst - integer(optional) Input: restart parameter +! + +subroutine psb_zgcr_vect(a,prec,b,x,eps,desc_a,info,& + & itmax,iter,err,itrace, irst, istop) + use psb_base_mod + use psb_prec_mod + use psb_z_krylov_conv_mod + use psb_krylov_mod + implicit none + + + type(psb_zspmat_type), intent(in) :: a + Type(psb_desc_type), Intent(in) :: desc_a + class(psb_zprec_type), intent(inout) :: prec + type(psb_z_vect_type), Intent(inout) :: b + type(psb_z_vect_type), Intent(inout) :: x + real(psb_dpk_), Intent(in) :: eps + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), Optional, Intent(in) :: itmax, itrace, irst, istop + integer(psb_ipk_), Optional, Intent(out) :: iter + real(psb_dpk_), Optional, Intent(out) :: err + ! = local data + complex(psb_dpk_), allocatable :: alpha(:), h(:,:) + type(psb_z_vect_type), allocatable :: z(:), c(:), c_scale(:) + type(psb_z_vect_type) :: r + + real(psb_dpk_) :: r_norm, b_norm, a_norm, derr + integer(psb_ipk_) :: n_col, mglob, naux, err_act + integer(psb_ipk_) :: debug_level, debug_unit + integer(psb_ipk_) :: np, me, ictxt + integer(psb_ipk_) :: i, j, it, itx, istop_, itmax_, itrace_, nl, m, nrst + complex(psb_dpk_) :: hjj + complex(psb_dpk_), allocatable, target :: aux(:) + character(len=20) :: name + type(psb_itconv_type) :: stopdat + character(len=*), parameter :: methdname='GCR' + integer(psb_ipk_) ::int_err(5) + info = psb_success_ + name = 'psb_zgcr' + call psb_erractionsave(err_act) + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + + ictxt = desc_a%get_context() + + call psb_info(ictxt, me, np) + if (.not.allocated(b%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + if (.not.allocated(x%v)) then + info = psb_err_invalid_vect_state_ + call psb_errpush(info,name) + goto 9999 + endif + + mglob = desc_a%get_global_rows() + n_col = desc_a%get_local_cols() + + if (present(istop)) then + istop_ = istop + else + istop_ = 2 + endif + + ! + ! ISTOP_ = 1: Normwise backward error, infinity norm + ! ISTOP_ = 2: ||r||/||b||, 2-norm + ! + + if ((istop_ < 1 ).or.(istop_ > 2 ) ) then + info=psb_err_invalid_istop_ + int_err(1)=istop_ + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + + call psb_chkvect(mglob,ione,x%get_nrows(),ione,ione,desc_a,info) + if (info == psb_success_)& + & call psb_chkvect(mglob,ione,b%get_nrows(),ione,ione,desc_a,info) + if(info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='psb_chkvect on X/B') + goto 9999 + end if + + if (present(itmax)) then + itmax_ = itmax + else + itmax_ = 1000 + endif + + if (present(itrace)) then + itrace_ = itrace + else + itrace_ = 0 + end if + + if (present(irst)) then + nl = irst + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' present: irst: ',irst,nl + else + nl = 10 + if (debug_level >= psb_debug_ext_) & + & write(debug_unit,*) me,' ',trim(name),& + & ' not present: irst: ',irst,nl + endif + + if (nl <=0 ) then + info=psb_err_invalid_istop_ + int_err(1)=nl + err=info + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif + + naux=4*n_col + allocate(aux(naux),h(nl+1,nl+1),& + &c_scale(nl+1),c(nl+1),z(nl+1), alpha(nl+1), stat=info) + + h = zzero + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_geasb(r, desc_a,info, scratch=.true.,mold=x%v) + + do i =1,nl+1 + call psb_geasb(c(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(z(i), desc_a,info, scratch=.true.,mold=x%v) + call psb_geasb(c_scale(i), desc_a,info, scratch=.true.,mold=x%v) + end do + + itx = 0 + + nrst = -1 + call psb_init_conv(methdname,istop_,itrace_,itmax_,a,x,b,eps,desc_a,stopdat,info) + restart: do + if (itx>= itmax_) exit restart + h = zzero + + it = 0 + ! compute r0 = b-ax0 + + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + + call psb_geaxpby(zone, b, zzero, r, desc_a, info) + call psb_spmm(-zone,a,x,zone,r,desc_a,info,work=aux) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + nrst = nrst + 1 + + iteration: do + + itx = itx + 1 + it = it + 1 + j = it + !Apply preconditioner + call prec%apply(r,z(j),desc_a,info,work=aux) + + call psb_spmm(zone,a,z(j),zzero,c(1),desc_a,info,work=aux) + do i =1, j - 1 + + h(i,j) = psb_gedot(c_scale(i), c(i), desc_a, info) + + call psb_geaxpby(zone, c(i), zzero, c(i+1), desc_a, info) + call psb_geaxpby(-h(i,j), c_scale(i), zone, c(i+1), desc_a, info) + end do + + h(j,j) = psb_norm2(c(j), desc_a, info) + hjj = zone/h(j,j) + call psb_geaxpby(hjj, c(j), zzero, c_scale(j), desc_a, info) + + alpha(j) = psb_gedot(c_scale(j), r, desc_a, info) + + !Update residual + call psb_geaxpby(zone, r, zzero, r, desc_a, info) + call psb_geaxpby(-alpha(j), c_scale(j), zone, r, desc_a, info) + + if (psb_check_conv(methdname,itx,x,r,desc_a,stopdat,info)) exit restart + + if (j >= irst) exit iteration + + + end do iteration + + m = j + + !Compute solution + + call ztrsm('l','u','n','n',m,1,zone,h,size(h,1),alpha,size(alpha,1)) + + if (nrst == 0 ) then + call x%set(zzero) + endif + do i=1,m + call psb_geaxpby(alpha(i), z(i), zone, x, desc_a, info) + enddo + + + + + end do restart + m = j + !Compute solution + call ztrsm('l','u','n','n',m,1,zone,h,size(h,1),alpha,size(alpha,1)) + call x%set(zzero) + do i=1,m + call psb_geaxpby(alpha(i), z(i), zone, x, desc_a, info) + enddo + + iter = j + + call psb_end_conv(methdname,itx,desc_a,stopdat,info,derr,iter) + if (present(err)) err = derr + + if (info == psb_success_) call psb_gefree(r,desc_a,info) + + do j = 1,m + if (info == psb_success_) call psb_gefree(z(j),desc_a,info) + if (info == psb_success_) call psb_gefree(c_scale(j),desc_a,info) + enddo + + do i =1,nl+1 + if (info == psb_success_) call psb_gefree(c(i),desc_a,info) + end do + + if (info == psb_success_) deallocate(aux,h,c_scale,z,c,alpha,stat=info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_non_ + call psb_errpush(info,name) + goto 9999 + end if + + call psb_erractionrestore(err_act) + return + +9999 continue + call psb_erractionrestore(err_act) + if (err_act.eq.psb_act_abort_) then + call psb_error() + return + end if + + return +end subroutine psb_zgcr_vect diff --git a/prec/CMakeLists.txt b/prec/CMakeLists.txt new file mode 100644 index 000000000..1fe2f5aa9 --- /dev/null +++ b/prec/CMakeLists.txt @@ -0,0 +1,66 @@ + +set(PSB_prec_source_files + psb_prec_const_mod.f90 + psb_c_prec_mod.f90 + psb_c_base_prec_mod.f90 + psb_d_bjacprec.f90 + psb_z_base_prec_mod.f90 + psb_d_diagprec.f90 + impl/psb_dprecbld.f90 + impl/psb_cprecinit.f90 + impl/psb_dprecinit.f90 + impl/psb_zilu_fct.f90 + impl/psb_silu_fct.f90 + impl/psb_cprecset.f90 + impl/psb_s_prec_type_impl.f90 + impl/psb_dilu_fct.f90 + impl/psb_cprecbld.f90 + impl/psb_c_diagprec_impl.f90 + impl/psb_s_bjacprec_impl.f90 + impl/psb_z_prec_type_impl.f90 + impl/psb_d_nullprec_impl.f90 + impl/psb_c_nullprec_impl.f90 + impl/psb_d_diagprec_impl.f90 + impl/psb_zprecinit.f90 + impl/psb_d_prec_type_impl.f90 + impl/psb_c_prec_type_impl.f90 + impl/psb_s_nullprec_impl.f90 + impl/psb_sprecinit.f90 + impl/psb_z_bjacprec_impl.f90 + impl/psb_d_bjacprec_impl.f90 + impl/psb_z_nullprec_impl.f90 + impl/psb_dprecset.f90 + impl/psb_c_bjacprec_impl.f90 + impl/psb_zprecset.f90 + impl/psb_sprecbld.f90 + impl/psb_cilu_fct.f90 + impl/psb_z_diagprec_impl.f90 + impl/psb_sprecset.f90 + impl/psb_zprecbld.f90 + impl/psb_s_diagprec_impl.f90 + psb_s_base_prec_mod.f90 + psb_s_bjacprec.f90 + psb_c_prec_type.f90 + psb_s_nullprec.f90 + psb_d_nullprec.f90 + psb_prec_mod.f90 + psb_d_prec_type.f90 + psb_z_prec_type.f90 + psb_z_prec_mod.f90 + psb_z_nullprec.f90 + psb_c_diagprec.f90 + psb_z_bjacprec.f90 + psb_s_prec_type.f90 + psb_c_nullprec.f90 + psb_prec_type.f90 + psb_d_prec_mod.f90 + psb_s_diagprec.f90 + psb_c_bjacprec.f90 + psb_s_prec_mod.f90 + psb_z_diagprec.f90 + psb_d_base_prec_mod.f90 + ) + +foreach(file IN LISTS PSB_prec_source_files) + list(APPEND prec_source_files ${CMAKE_CURRENT_LIST_DIR}/${file}) +endforeach() diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt new file mode 100644 index 000000000..93464741a --- /dev/null +++ b/test/CMakeLists.txt @@ -0,0 +1,94 @@ +if (NOT MPI_C_FOUND) + find_package(MPI REQUIRED) + + set(CMAKE_C_COMPILE_FLAGS ${CMAKE_C_COMPILE_FLAGS} ${MPI_C_COMPILE_FLAGS}) + set(CMAKE_C_LINK_FLAGS ${CMAKE_C_LINK_FLAGS} ${MPI_C_LINK_FLAGS}) + set(CMAKE_Fortran_COMPILE_FLAGS ${CMAKE_Fortran_COMPILE_FLAGS} ${MPI_Fortran_COMPILE_FLAGS}) + set(CMAKE_Fortran_LINK_FLAGS ${CMAKE_Fortran_LINK_FLAGS} ${MPI_Fortran_LINK_FLAGS}) + include_directories(BEFORE ${MPI_C_INCLUDE_PATH} ${MPI_Fortran_INCLUDE_PATH}) +endif() + +if("${CMAKE_Fortran_COMPILER_ID}" STREQUAL "GNU") + set(gfortran_compiler true) +endif() + +include(CheckIncludeFile) +CHECK_INCLUDE_FILE("alloca.h" HAVE_ALLOCA) +if(NOT HAVE_ALLOCA) + add_definitions(-DALLOCA_MISSING) + message(WARNING "Could not find . Assuming functionality is provided elsewhere.") +endif() + +set(old_cmake_required_includes "${CMAKE_REQUIRED_INCLUDES}") +if(CMAKE_REQUIRED_INCLUDES) + set(CMAKE_REQUIRED_INCLUDES ${CMAKE_REQUIRED_INCLUDES};${MPI_C_INCLUDE_PATH}) +else() + set(CMAKE_REQUIRED_INCLUDES ${MPI_C_INCLUDE_PATH}) +endif() +set(old_cmake_required_flags "${CMAKE_REQUIRED_FLAGS}") +if(CMAKE_REQUIRED_FLAGS) + set(CMAKE_REQUIRED_FLAGS ${CMAKE_REQUIRED_FLAGS};${MPI_C_COMPILE_FLAGS};${MPI_C_LINK_FLAGS}) +else() + set(CMAKE_REQUIRED_FLAGS ${MPI_C_COMPILE_FLAGS};${MPI_C_LINK_FLAGS}) +endif() +set(old_cmake_required_libraries "${CMAKE_REQUIRED_LIBRARIES}") +if(CMAKE_REQUIRED_LIBRARIES) + set(CMAKE_REQUIRED_LIBRARIES ${CMAKE_REQUIRED_LIBRARIES};${MPI_C_LIBRARIES}) +else() + set(CMAKE_REQUIRED_LIBRARIES ${MPI_C_LIBRARIES}) +endif() + +set(MPI_HEADERS mpi.h) +include(CheckIncludeFiles) +CHECK_INCLUDE_FILES("mpi.h;mpi-ext.h" HAVE_MPI_EXT) +if(HAVE_MPI_EXT) + add_definitions(-DHAVE_MPI_EXT_H) + set(MPI_HEADERS ${MPI_HEADERS};mpi-ext.h) +endif() + +set(CMAKE_REQUIRED_INCLUDES ${old_cmake_required_includes}) +set(CMAKE_REQUIRED_FLAGS ${old_cmake_required_flags}) +set(CMAKE_REQUIRED_LIBRARIES ${old_cmake_required_libraries}) + +#--------------------------------------------------- +# Windows Intel MPI compatibility, see GH issue #435 +#--------------------------------------------------- +set(old_cmake_required_includes "${CMAKE_REQUIRED_INCLUDES}") +if(CMAKE_REQUIRED_INCLUDES) + set(CMAKE_REQUIRED_INCLUDES ${CMAKE_REQUIRED_INCLUDES};${MPI_C_INCLUDE_PATH}) +else() + set(CMAKE_REQUIRED_INCLUDES ${MPI_C_INCLUDE_PATH}) +endif() +CHECK_INCLUDE_FILES("mpi.h" HAVE_MPI_H) + +include(CheckSymbolExists) +CHECK_SYMBOL_EXISTS(I_MPI_VERSION "mpi.h" HAVE_Intel_MPI) +if(HAVE_Intel_MPI AND WIN32) + add_definitions(-DUSE_GCC) +endif() +set(CMAKE_REQUIRED_INCLUDES ${old_cmake_required_includes}) + +set(PSBLAS_SO_VERSION 0) +if(gfortran_compiler) + if(NOT CMAKE_Fortran_COMPILER_VERSION VERSION_LESS 7.0.0) + set(PSBLAS_SO_VERSION 2) + elseif(NOT CMAKE_Fortran_COMPILER_VERSION VERSION_LESS 6.0.0) + set(PSBLAS_SO_VERSION 1) + endif() +endif() + +set(PSBLAS_MODDIR "${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") +set(MOD_DIR_FLAG "${CMAKE_Fortran_MODDIR_FLAG}") +set(MPI_LIBS "") +foreach( lib IN LISTS MPI_Fortran_LIBRARIES) + set(MPI_LIBS "${MPI_LIBS} \"${lib}\"") +endforeach() +string(STRIP "${MPI_LIBS}" MPI_LIBS) + +set(directory_list hello )#fileread hello kernel pargen torture util) +if (SERIAL_MPI) + set(directory_list ${directory_list} serial) +endif() +foreach(directory ${directory_list}) + add_subdirectory("${directory}") +endforeach() diff --git a/test/fileread/CMakeLists.txt b/test/fileread/CMakeLists.txt new file mode 100644 index 000000000..4362687d0 --- /dev/null +++ b/test/fileread/CMakeLists.txt @@ -0,0 +1,11 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable(psb_cf_sample psb_cf_sample.f90 getp.f90) +add_executable(psb_df_sample psb_df_sample.f90 getp.f90) +add_executable(psb_sf_sample psb_sf_sample.f90 getp.f90) +add_executable(psb_zf_sample psb_zf_sample.f90 getp.f90) + +target_link_libraries(psb_cf_sample krylov util base) +target_link_libraries(psb_df_sample krylov util base) +target_link_libraries(psb_sf_sample krylov util base) +target_link_libraries(psb_zf_sample krylov util base) diff --git a/test/hello/CMakeLists.txt b/test/hello/CMakeLists.txt new file mode 100644 index 000000000..e3884fad3 --- /dev/null +++ b/test/hello/CMakeLists.txt @@ -0,0 +1,7 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable(hello hello.f90) +add_executable(pingpong pingpong.f90) + +target_link_libraries(hello base) +target_link_libraries(pingpong base) diff --git a/test/kernel/CMakeLists.txt b/test/kernel/CMakeLists.txt new file mode 100644 index 000000000..f80eff8b7 --- /dev/null +++ b/test/kernel/CMakeLists.txt @@ -0,0 +1,9 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable(d_file_spmv d_file_spmv.f90) +add_executable(pdgenspmv pdgenspmv.f90) +add_executable(s_file_spmv s_file_spmv.f90) + +target_link_libraries(d_file_spmv util) +target_link_libraries(pdgenspmv util) +target_link_libraries(s_file_spmv util) diff --git a/test/pargen/CMakeLists.txt b/test/pargen/CMakeLists.txt new file mode 100644 index 000000000..443c1f6ef --- /dev/null +++ b/test/pargen/CMakeLists.txt @@ -0,0 +1,11 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable(psb_d_pde2d psb_d_pde2d.f90) +add_executable(psb_d_pde3d psb_d_pde3d.f90) +add_executable(psb_s_pde2d psb_s_pde2d.f90) +add_executable(psb_s_pde3d psb_s_pde3d.f90) + +target_link_libraries(psb_d_pde2d krylov util) +target_link_libraries(psb_d_pde3d krylov util) +target_link_libraries(psb_s_pde2d krylov util) +target_link_libraries(psb_s_pde3d krylov util) diff --git a/test/serial/CMakeLists.txt b/test/serial/CMakeLists.txt new file mode 100644 index 000000000..3556f8b39 --- /dev/null +++ b/test/serial/CMakeLists.txt @@ -0,0 +1,8 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable( d_matgen + psb_d_xyz_impl.f90 + psb_d_xyz_mat_mod.f90 + d_matgen.F90 # Main program +) +target_link_libraries(d_matgen base) diff --git a/test/torture/CMakeLists.txt b/test/torture/CMakeLists.txt new file mode 100644 index 000000000..7c5c27517 --- /dev/null +++ b/test/torture/CMakeLists.txt @@ -0,0 +1,11 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable(fortran_interface_test + psb_c_mvsv_tester.f90 + psb_d_mvsv_tester.f90 + psb_mvsv_tester.f90 + psb_s_mvsv_tester.f90 + psbtf.f90 # Main program + psb_z_mvsv_tester.f90 +) +target_link_libraries(fortran_interface_test base) diff --git a/test/util/CMakeLists.txt b/test/util/CMakeLists.txt new file mode 100644 index 000000000..a38e54c2a --- /dev/null +++ b/test/util/CMakeLists.txt @@ -0,0 +1,11 @@ +include_directories("${CMAKE_BINARY_DIR}/${CMAKE_INSTALL_INCLUDEDIR}/${mod_dir_tail}") + +add_executable(dhb2mm dhb2mm.f90) +add_executable(dmm2hb dmm2hb.f90) +add_executable(zhb2mm zhb2mm.f90) +add_executable(zmm2hb zmm2hb.f90) + +target_link_libraries(dhb2mm util) +target_link_libraries(dmm2hb util) +target_link_libraries(zhb2mm util) +target_link_libraries(zmm2hb util) diff --git a/test/util/dhb2mm.f90 b/test/util/dhb2mm.f90 index 1e66608b1..e31b143c1 100644 --- a/test/util/dhb2mm.f90 +++ b/test/util/dhb2mm.f90 @@ -1,9 +1,9 @@ -! +! ! Parallel Sparse BLAS version 3.5 ! (C) Copyright 2006, 2010, 2015, 2017 ! Salvatore Filippone Cranfield University ! Alfredo Buttari CNRS-IRIT, Toulouse -! +! ! Redistribution and use in source and binary forms, with or without ! modification, are permitted provided that the following conditions ! are met: @@ -15,7 +15,7 @@ ! 3. The name of the PSBLAS group or the names of its contributors may ! not be used to endorse or promote products derived from this ! software without specific 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 @@ -27,22 +27,22 @@ ! 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. -! -! ! -! Storage conversion filter: reads from standar input a sparse matrix -! stored in Harwell-Boeing format, and writes to standard output in MatrixMarket +! +! +! Storage conversion filter: reads from standard input a sparse matrix +! stored in Harwell-Boeing format, and writes to standard output in MatrixMarket ! format ! program dhb2mm use psb_base_mod use psb_util_mod type(psb_dspmat_type) :: a - + integer(psb_ipk_) :: info character(len=72) :: mtitle - + call hb_read(a,info,mtitle=mtitle) call mm_mat_write(a,mtitle,info) @@ -51,4 +51,4 @@ program dhb2mm end program dhb2mm - + diff --git a/util/CMakeLists.txt b/util/CMakeLists.txt new file mode 100644 index 000000000..5c77191e7 --- /dev/null +++ b/util/CMakeLists.txt @@ -0,0 +1,42 @@ +set(PSB_util_source_files + psb_blockpart_mod.f90 + psb_c_hbio_impl.f90 + psb_c_mat_dist_impl.f90 + psb_c_mat_dist_mod.f90 + psb_c_mmio_impl.f90 + psb_c_renum_impl.F90 + psb_d_hbio_impl.f90 + psb_d_mat_dist_impl.f90 + psb_d_mat_dist_mod.f90 + psb_d_mmio_impl.f90 + psb_d_renum_impl.F90 + psb_gps_mod.f90 + psb_hbio_mod.f90 + psb_i_mmio_impl.f90 + psb_mat_dist_mod.f90 + psb_metispart_mod.F90 + psb_mmio_mod.F90 + psb_renum_mod.f90 + psb_s_hbio_impl.f90 + psb_s_mat_dist_impl.f90 + psb_s_mat_dist_mod.f90 + psb_s_mmio_impl.f90 + psb_s_renum_impl.F90 + psb_util_mod.f90 + psb_z_hbio_impl.f90 + psb_z_mat_dist_impl.f90 + psb_z_mat_dist_mod.f90 + psb_z_mmio_impl.f90 + psb_z_renum_impl.F90 + ) +foreach(file IN LISTS PSB_util_source_files) + list(APPEND util_source_files ${CMAKE_CURRENT_LIST_DIR}/${file}) +endforeach() + +set(PSB_util_source_C_files + metis_int.c + psb_amd_order.c + ) +foreach(file IN LISTS PSB_util_source_C_files) + list(APPEND util_source_C_files ${CMAKE_CURRENT_LIST_DIR}/${file}) +endforeach()