diff --git a/fem/src/DirectSolve.F90 b/fem/src/DirectSolve.F90 index 2e0916ce6bf8384e237ad836aab4c2334cb6d9e9..a7f13954325218858dbecc868fe31a603da865db 100644 --- a/fem/src/DirectSolve.F90 +++ b/fem/src/DirectSolve.F90 @@ -1883,8 +1883,8 @@ CONTAINS INTEGER :: n_dof_partition - INTEGER, ALLOCATABLE :: memb(:), DirichletInds(:), Neighbours(:) - REAL(KIND=dp), ALLOCATABLE :: DirichletVals(:) + INTEGER, ALLOCATABLE :: memb(:), DirichletInds(:), Neighbours(:), IntOptions(:) + REAL(KIND=dp), ALLOCATABLE :: DirichletVals(:), RealOptions(:) INTEGER :: Comm_active, Group_active, Group_world REAL(KIND=dp), ALLOCATABLE :: dbuf(:) @@ -1950,10 +1950,17 @@ CONTAINS 'inconsistency: A % NumberOfRows /= SIZE(A % ParallelInfo % GlobalDOFs' ) END IF + ALLOCATE(IntOptions(10)) + ALLOCATE(RealOptions(1)) + + CALL FETI4ISetDefaultIntegerOptions(IntOptions) + CALL FETI4ISetDefaultRealOptions(RealOptions) + CALL FETI4ICreateInstance(A % PermonSolverInstance, A % PermonMatrix, & A % NumberOfRows, b, A % ParallelInfo % GlobalDOFs, & n, neighbours, & - nd, DirichletInds, DirichletVals) + nd, DirichletInds, DirichletVals, & + IntOptions, RealOptions) END IF !CALL Permon_Solve( A % PermonSolverInstance, x, b ) diff --git a/feti4i/feti4i.h b/feti4i/feti4i.h deleted file mode 100644 index 1423e6ce1b9ba738bafd4e3dc80c4522d2f1e103..0000000000000000000000000000000000000000 --- a/feti4i/feti4i.h +++ /dev/null @@ -1,134 +0,0 @@ - -#ifndef FETI4I_H_ -#define FETI4I_H_ - -#include "mpi.h" - -/*----------------------------------------------------------------------------- - Set data-types used in FETI4I - - Possible values for FETI4I_INT_WIDTH: - 32: use 32 bit signed integer - 64: use 64 bit signed integer - - Possible values for FETI4I_REAL_WIDTH: - 64: ESPRESO supports only 64 bit real values -------------------------------------------------------------------------------*/ -#ifndef FETI4I_INT_WIDTH -#define FETI4I_INT_WIDTH 32 -#endif - -#ifndef FETI4I_REAL_WIDTH -#define FETI4I_REAL_WIDTH 64 -#endif - -#if FETI4I_INT_WIDTH == 32 - typedef int FETI4IInt; -#elif FETI4I_INT_WIDTH == 64 - typedef long FETI4IInt; -#else -#error "Incorrect user-supplied value of FETI4I_INT_WIDTH" -#endif - -/* MPI integer (e.g. rank) are always 32-bit */ - typedef int FETI4IMPIInt; - -#if FETI4I_REAL_WIDTH == 64 - typedef double FETI4IReal; -#else -#error "Incorrect user-supplied value of FETI4I_REAL_WIDTH" -#endif - - -#ifdef __cplusplus -extern "C" { -#endif - -/*----------------------------------------------------------------------------- - Definitions of internal structures used in FETI4I -------------------------------------------------------------------------------*/ -typedef struct FETI4IStructMatrix* FETI4IMatrix; -typedef struct FETI4IStructInstance* FETI4IInstance; - -/*----------------------------------------------------------------------------- - Functions for manipulating with FETI4I internal structures -------------------------------------------------------------------------------*/ - -void FETI4ICreateStiffnessMatrix( - FETI4IMatrix *matrix, //TODO size? - FETI4IInt indexBase -); - -void FETI4IAddElement( //TODO CRS or CCS ? - FETI4IMatrix matrix, - FETI4IInt size, - FETI4IInt* indices, - FETI4IReal* values -); - -/*----------------------------------------------------------------------------- - Functions for creating an instance and solve it -------------------------------------------------------------------------------*/ - -void FETI4ICreateInstance( - FETI4IInstance *instance, - FETI4IMatrix matrix, - FETI4IInt size, - FETI4IReal* rhs, - FETI4IInt* l2g, /* length of both rhs and l2g is size */ - FETI4IMPIInt neighbours_size, - FETI4IMPIInt* neighbours, - FETI4IInt dirichlet_size, - FETI4IInt* dirichlet_indices, //TODO which numbering? we prefer global numbering - FETI4IReal* dirichlet_values -); - -void FETI4ISolve( - FETI4IInstance instance, - FETI4IInt solution_size, - FETI4IReal* solution -); - -/*----------------------------------------------------------------------------- - Functions for updating a created instance -------------------------------------------------------------------------------*/ - -void FETI4IUpdateStiffnessMatrix( - FETI4IInstance instance, - FETI4IMatrix stiffnessMatrix -); - -void FETI4IUpdateRhs( - FETI4IInstance instance, - FETI4IInt rhs_size, - FETI4IReal* rhs_values -); - -//TODO VH: neighbours perhaps should not be passed here and they are FETI4IMPIInt -void FETI4IUpdateDirichlet( - FETI4IInstance instance, - FETI4IInt dirichlet_size, - FETI4IInt* dirichlet_indices, - FETI4IReal* dirichlet_values, - FETI4IReal* l2g, - FETI4IInt neighbour_size, - FETI4IReal* neighbour -); - - -/*----------------------------------------------------------------------------- - Destroy an arbitrary internal structure -------------------------------------------------------------------------------*/ - -void FETI4IDestroy( - void* ptr -); - - -#ifdef __cplusplus -} -#endif - - -#endif /* FETI4I_H_ */ - diff --git a/feti4i/feti4i_fortran_def.h90 b/feti4i/feti4i_fortran_def.h90 deleted file mode 100644 index a688bab4ab2abd20d1b96f58174bd604c4de9bf6..0000000000000000000000000000000000000000 --- a/feti4i/feti4i_fortran_def.h90 +++ /dev/null @@ -1,50 +0,0 @@ -!#ifndef FETI4I_FORTRAN_DEF_H_ -!#define FETI4I_FORTRAN_DEF_H_ - -! ------------------------------------------------------------------------------ -! Set data-types used in FETI4I -! -! Possible values for FETI4I_INT_WIDTH: -! 32: use 32 bit signed integer -! 64: use 64 bit signed integer -! -! Possible values for FETI4I_REAL_WIDTH: -! 64: ESPRESO supports only 64 bit real values -! ------------------------------------------------------------------------------ - -USE, INTRINSIC :: ISO_C_BINDING - - -#ifndef FETI4I_INT_WIDTH -#define FETI4I_INT_WIDTH 32 -#endif - -#ifndef FETI4I_REAL_WIDTH -#define FETI4I_REAL_WIDTH 64 -#endif - -#if FETI4I_INT_WIDTH == 32 - INTEGER, PARAMETER :: FETI4IInt = C_INT -#elif FETI4I_INT_WIDTH == 64 - INTEGER, PARAMETER :: FETI4IInt = C_LONG -#else -#error "Incorrect user-supplied value of FETI4I_INT_WIDTH" -#endif - -! MPI integer (e.g. rank) are always 32-bit - INTEGER, PARAMETER :: FETI4IMPIInt = C_INT - -#if FETI4I_REAL_WIDTH == 64 - INTEGER, PARAMETER :: FETI4IReal = C_DOUBLE -#else -#error "Incorrect user-supplied value of FETI4I_REAL_WIDTH" -#endif - -! ----------------------------------------------------------------------------- -! Definitions of internal structures used in FETI4I -! ----------------------------------------------------------------------------- -! use TYPE(C_PTR) for opaque pointers -!TODO we could perhaps use Fortran types - -!#endif /* FETI4I_FORTRAN_DEF_H_ */ - diff --git a/feti4i/feti4i_fortran_interface.h90 b/feti4i/feti4i_fortran_interface.h90 deleted file mode 100644 index fbb97da1c3449bd08ff1e25c78e8b33866e0bba5..0000000000000000000000000000000000000000 --- a/feti4i/feti4i_fortran_interface.h90 +++ /dev/null @@ -1,52 +0,0 @@ - -#include "feti4i_fortran_def.h90" - -INTERFACE - SUBROUTINE FETI4ICreateStiffnessMatrix(matrix, indexBase) & - BIND(c,NAME='FETI4ICreateStiffnessMatrix') - IMPORT - TYPE(C_PTR) :: matrix ! FETI4IMatrix* type in C (output) - INTEGER(KIND=FETI4IInt), VALUE :: indexBase - END SUBROUTINE FETI4ICreateStiffnessMatrix - - SUBROUTINE FETI4IAddElement(matrix, size, indices, values) & - BIND(c,NAME='FETI4IAddElement') - IMPORT - TYPE(C_PTR), VALUE :: matrix ! FETI4IMatrix type in C (input) - INTEGER(KIND=FETI4IInt), VALUE :: size - INTEGER(KIND=FETI4IInt) :: indices(*) - REAL(KIND=FETI4IReal) :: values(*) - END SUBROUTINE FETI4IAddElement - - SUBROUTINE FETI4ICreateInstance(instance, matrix, size, rhs, l2g, & - neighbours_size, neighbours, & - dirichlet_size, dirichlet_indices, dirichlet_values) & - BIND(c,NAME='FETI4ICreateInstance') - IMPORT - TYPE(C_PTR) :: instance ! FETI4IInstance* type in C (output) - TYPE(C_PTR), VALUE :: matrix ! FETI4IMatrix type in C (input) - INTEGER(KIND=FETI4IInt), VALUE :: size - REAL(KIND=FETI4IReal) :: rhs(*) - INTEGER(KIND=FETI4IInt) :: l2g(*) - INTEGER(KIND=FETI4IMPIInt), VALUE :: neighbours_size - INTEGER(KIND=FETI4IMPIInt) :: neighbours(*) - INTEGER(KIND=FETI4IInt), VALUE :: dirichlet_size - INTEGER(KIND=FETI4IInt) :: dirichlet_indices(*) - REAL(KIND=FETI4IReal) :: dirichlet_values(*) - END SUBROUTINE FETI4ICreateInstance - - SUBROUTINE FETI4ISolve(instance, solution_size, solution) & - BIND(c,NAME='FETI4ISolve') - IMPORT - TYPE(C_PTR), VALUE :: instance ! FETI4IInstance type in C (input) - INTEGER(KIND=FETI4IInt), VALUE :: solution_size - REAL(KIND=FETI4IReal) :: solution(*) - END SUBROUTINE FETI4ISolve - - SUBROUTINE FETI4IDestroy(ptr) & - BIND(c,NAME='FETI4IDestroy') - IMPORT - TYPE(C_PTR), VALUE :: ptr ! void* type in C - END SUBROUTINE FETI4IDestroy -END INTERFACE - diff --git a/feti4i/feti4i_fortran_test/Makefile b/feti4i/feti4i_fortran_test/Makefile deleted file mode 100644 index 6e15fa43c2eecf6cbafe0ee7021992d36ac92568..0000000000000000000000000000000000000000 --- a/feti4i/feti4i_fortran_test/Makefile +++ /dev/null @@ -1,35 +0,0 @@ -FC := mpif90 -CC := mpicc -MKDIR := mkdir -p -OUTLIBDIR := ../lib - -build: test - -feti4i_mod.o: ../feti4i_mod.F90 - ${FC} -c -fPIC $^ - -feti4i.mod: ../feti4i_mod.F90 feti4i_mod.o - @true - -feti4i_dummy.o: feti4i_dummy.c - ${CC} -c -fPIC -g $^ - -${OUTLIBDIR}/libfeti4i.so: feti4i_dummy.o - ${MKDIR} ${OUTLIBDIR} - ${LINK.F} -o ${OUTLIBDIR}/libfeti4i.so feti4i_dummy.o -shared -fPIC - -test.o: test.F90 feti4i.mod - ${FC} -c -fPIC -g test.F90 -o $@ - -test: test.o ${OUTLIBDIR}/libfeti4i.so - ${LINK.F} -o test test.o -L${OUTLIBDIR} -Wl,-rpath=${OUTLIBDIR} -lfeti4i - -clean: - ${RM} test.o feti4i_dummy.o feti4i_mod.o - -clean-all: clean - ${RM} feti4i.mod test - ${RM} -r ${OUTLIBDIR} - -.PHONY: build clean - diff --git a/feti4i/feti4i_fortran_test/feti4i_dummy.c b/feti4i/feti4i_fortran_test/feti4i_dummy.c deleted file mode 100644 index 7e160afeac4e3a5a9d7e981e9716c2def0c5bdef..0000000000000000000000000000000000000000 --- a/feti4i/feti4i_fortran_test/feti4i_dummy.c +++ /dev/null @@ -1,95 +0,0 @@ - -#include <stdlib.h> -#include <stdio.h> -#include <assert.h> -#include "../feti4i.h" - -struct FETI4IStructMatrix { - FETI4IInt indexBase; - FETI4IInt size; - FETI4IInt* indices; - FETI4IReal* values; - FETI4IInt nelems; -}; - -struct FETI4IStructInstance { - FETI4IMatrix matrix; - FETI4IReal* rhs; - FETI4IInt* l2g; /* length of both rhs and l2g is size */ - FETI4IMPIInt neighbours_size; - FETI4IMPIInt* neighbours; - FETI4IInt dirichlet_size; - FETI4IInt* dirichlet_indices; //TODO which numbering? we prefer global numbering - FETI4IReal* dirichlet_values; -}; - - -void FETI4ICreateStiffnessMatrix(FETI4IMatrix *matrix, FETI4IInt indexBase) -{ - *matrix = (FETI4IMatrix) malloc(sizeof(struct FETI4IStructMatrix)); - (*matrix)->indexBase = indexBase; - (*matrix)->size = 1; - (*matrix)->indices = NULL; - (*matrix)->values = NULL; - (*matrix)->nelems = 0; - printf("%s: created matrix %p, indexBase %d\n", __func__, (void*) *matrix, indexBase); -} - - -void FETI4IAddElement(FETI4IMatrix matrix, FETI4IInt size, FETI4IInt* indices, FETI4IReal* values) -{ - matrix->size += (size-1); - matrix->indices = indices; - matrix->values = values; - printf("%s: added element %d to matrix %p\n", __func__, matrix->nelems++, (void*) matrix); -} - -void FETI4ICreateInstance(FETI4IInstance *instance, FETI4IMatrix matrix, FETI4IInt size, - FETI4IReal* rhs, FETI4IInt* l2g, FETI4IMPIInt neighbours_size, FETI4IMPIInt* neighbours, - FETI4IInt dirichlet_size, FETI4IInt* dirichlet_indices, FETI4IReal* dirichlet_values) -{ - printf("%s: setting rhs %p (%d), l2g %p (%d), " - "neighbours %p (%d), dirichlet_indices %p (%d), " - "dirichlet_values %p (%d)\n", - __func__, (void*) rhs, size, (void*) l2g, size, - (void*) neighbours, neighbours_size, (void*) dirichlet_indices, dirichlet_size, - (void*) dirichlet_indices, dirichlet_size); - - *instance = (FETI4IInstance) malloc(sizeof(struct FETI4IStructInstance)); - (*instance)->matrix = matrix; - //assert(size == matrix->size); - (*instance)->rhs = rhs; - (*instance)->l2g = l2g; - (*instance)->neighbours_size = neighbours_size; - (*instance)->neighbours = neighbours; - (*instance)->dirichlet_size = dirichlet_size; - (*instance)->dirichlet_indices = dirichlet_indices; - (*instance)->dirichlet_values = dirichlet_values; -} - -void FETI4ISolve(FETI4IInstance instance, FETI4IInt solution_size, FETI4IReal* solution) -{ - FETI4IInt i; - FETI4IMatrix matrix = instance->matrix; - FETI4IReal *rhs = instance->rhs; - - printf("%s: solving with matrix %p, rhs %p\n", __func__, (void*) instance->matrix, (void*) instance->rhs); - //assert(solution_size == matrix->size); - //for (i=0; i<solution_size; i++) { - // printf("#%d %g\n",i,rhs[i]); - //} - //printf("\n"); - //for (i=0; i<solution_size; i++) { - // printf("#%d %g ",i,solution[i]); - // solution[i]=-i; - // printf(" %g\n",solution[i]); - //} -} - -void FETI4IDestroy(void *ptr) -{ - FETI4IInstance instance = (FETI4IInstance) ptr; - free(instance->matrix); - free(instance); -} - diff --git a/feti4i/feti4i_fortran_test/test.F90 b/feti4i/feti4i_fortran_test/test.F90 deleted file mode 100644 index c3beb7e8e5c39ca849d712ca0bed5c4261a5273d..0000000000000000000000000000000000000000 --- a/feti4i/feti4i_fortran_test/test.F90 +++ /dev/null @@ -1,74 +0,0 @@ -PROGRAM test - - USE FETI4I - - TYPE(C_PTR) :: instance - TYPE(C_PTR) :: matrix - - INTEGER(KIND=FETI4IInt) :: n = 4 - INTEGER(KIND=FETI4IInt) :: nelems = 5 - INTEGER(KIND=FETI4IInt) :: size = 0 !TODO when this is known - INTEGER(KIND=FETI4IInt), ALLOCATABLE :: inds(:) - REAL(KIND=FETI4IReal), ALLOCATABLE :: vals(:) - REAL(KIND=FETI4IReal), ALLOCATABLE :: b(:) - INTEGER(KIND=FETI4IInt), ALLOCATABLE :: l2g(:) - INTEGER(KIND=FETI4IMPIInt) :: neighbours_size = 3 - INTEGER(KIND=FETI4IMPIInt), ALLOCATABLE :: neighbours(:) - INTEGER(KIND=FETI4IInt) :: dirichlet_size = 2 - INTEGER(KIND=FETI4IInt), ALLOCATABLE :: dirichlet_indices(:) - REAL(KIND=FETI4IReal), ALLOCATABLE :: dirichlet_values(:) - REAL(KIND=FETI4IReal), ALLOCATABLE :: x(:) - - INTEGER :: ALLOC_ERR - - CALL FETI4ICreateStiffnessMatrix(matrix, 1) - - ALLOCATE(vals(n*n), inds(n)) - DO i=1,n - DO j=1,n - vals((i-1)*n+j) = (i-1)*n+j - END DO - inds(i) = i - END DO - - size = 1 - DO i=1,nelems - CALL FETI4IAddElement(matrix, n, inds, vals) - size = size + (n-1) - END DO - - ALLOCATE(b(size), l2g(size)) - DO i=1,size - b(i) = -i - l2g(i) = 100+i - END DO - - write (*,*) "b: ", b - - ALLOCATE(neighbours(neighbours_size)) - DO i=1,neighbours_size - neighbours(i) = i-1 - END DO - - ALLOCATE(dirichlet_indices(dirichlet_size), dirichlet_values(dirichlet_size)) - DO i=1,dirichlet_size - dirichlet_indices(i) = 2*(i+1) - dirichlet_values(i) = -2*(i+1) - END DO - - CALL FETI4ICreateInstance(instance, matrix, size, b, l2g, & - neighbours_size, neighbours, & - dirichlet_size, dirichlet_indices, dirichlet_values) - - ALLOCATE(x(size), STAT = ALLOC_ERR) - write (*,*) "ALLOC_ERR=",ALLOC_ERR - - x(:) = 1.0 - - write (*,*) "x: ", x - - CALL FETI4ISolve(instance, size, x) - - CALL FETI4IDestroy(instance) - -END PROGRAM test diff --git a/feti4i/feti4i_mod.F90 b/feti4i/feti4i_mod.F90 deleted file mode 100644 index ac14564402f43135649f52d5f8c58361aa514553..0000000000000000000000000000000000000000 --- a/feti4i/feti4i_mod.F90 +++ /dev/null @@ -1,5 +0,0 @@ - -MODULE feti4i -#include "feti4i_fortran_interface.h90" -END MODULE feti4i -