!DEC$ FREEFORM
!=============================================================================
! Monolithic FE^2
! Nils Lange, Geralf Huetter, Bjoern Kiefer
!   Nils.Lange@imfd.tu-freiberg.de, Geralf.Huetter@imfd.tu-freiberg.de, 
!   Bjoern.Kiefer@imfd.tu-freiberg.de
! distributed under CC BY-NC-SA 4.0 license
! (https://creativecommons.org/licenses/by-nc-sa/4.0/)
! Reference: 
!   N. Lange, G. Huetter, B. Kiefer: "An efficient monolithic solution scheme for FE2 problems",
!   https://arxiv.org/abs/2101.01802
!
! Further information on the implementation, structure of the source code,
! examples and tutorials can be found in the file doc/documentation.pdf
! 
!=============================================================================

INCLUDE 'Solver_LAPACK.f'                 !LAPACK Solver Module
INCLUDE 'Solver_PARDISO.f'                !PARDISO Solver Module

MODULE class_Solver

!general "empty" interface for any solver; which actual solver will be used is
!determined by the value from the analysisparameters value of solving_process

    USE Solver_LAPACK
    USE Solver_PARDISO
    
    TYPE solve_soe
    
        TYPE(solve_soe_PARDISO):: PARDISOsolver
        TYPE(solve_soe_LAPACK)::  LAPACKsolver
        
        CONTAINS
        
            PROCEDURE :: initialize => initialize_procedure
            PROCEDURE :: factor => factor_procedure
            PROCEDURE :: solve => solve_procedure
            PROCEDURE :: finish => finish_procedure
    
    END TYPE solve_soe
    
    PRIVATE:: initialize_procedure,factor_procedure,solve_procedure,finish_procedure
    
    CONTAINS
    
        SUBROUTINE initialize_procedure(solver,para,analysis,GP_DATA)
        
            IMPLICIT NONE
        
            CLASS(solve_soe), INTENT(INOUT):: solver
            TYPE (analysisparameters)::analysis
            TYPE(meshparameters),INTENT(IN)::para
            TYPE(macro_GP_DATA),INTENT(INOUT)::GP_DATA
            
            IF (analysis%solving_process==0) THEN !full simulation
                CALL solver%PARDISOsolver%initialize(para,analysis,GP_DATA)
            ELSE IF (analysis%solving_process==1) THEN !reduced simulation
                CALL solver%LAPACKsolver%initialize(para,analysis,GP_DATA)
            END IF
                
        END SUBROUTINE initialize_procedure
        
        SUBROUTINE factor_procedure(solver,para,analysis,k_matrix_values,k_matrix)
        
            IMPLICIT NONE
        
            CLASS(solve_soe),INTENT(INOUT):: solver
            TYPE (analysisparameters),INTENT(IN)::analysis
            TYPE(meshparameters),INTENT(IN)::para
            REAL(KIND=AbqRK),DIMENSION(:),INTENT(IN)::k_matrix_values
            REAL(KIND=AbqRK),DIMENSION(:,:),INTENT(IN),OPTIONAL::k_matrix
            
            IF (analysis%solving_process==0) THEN !full simulation
                CALL solver%PARDISOsolver%factor(para,analysis,k_matrix_values)
            ELSE IF (analysis%solving_process==1) THEN !reduced simulation
                CALL solver%LAPACKsolver%factor(para,analysis,k_matrix)
            END IF
                
        END SUBROUTINE factor_procedure
        
        SUBROUTINE solve_procedure(solver,para,analysis,k_matrix_values,rhs,solution)
        
            IMPLICIT NONE
        
            CLASS(solve_soe),INTENT(INOUT):: solver
            TYPE (analysisparameters),INTENT(IN)::analysis
            TYPE(meshparameters),INTENT(IN)::para
            REAL(KIND=AbqRK),DIMENSION(:),INTENT(IN)::k_matrix_values
            REAL(KIND=AbqRK),DIMENSION(:,:),INTENT(IN):: rhs
            REAL(KIND=AbqRK),DIMENSION(:,:),INTENT(OUT):: solution
            
            IF (analysis%solving_process==0) THEN !full simulation
                CALL solver%PARDISOsolver%solve(para,analysis,k_matrix_values,rhs,solution)
            ELSE IF (analysis%solving_process==1) THEN !reduced simulation
                CALL solver%LAPACKsolver%solve(para,analysis,rhs,solution)
            END IF
                
        END SUBROUTINE solve_procedure
    
        SUBROUTINE finish_procedure(solver,para,analysis,GP_DATA)
        
            IMPLICIT NONE
        
            CLASS(solve_soe),INTENT(INOUT):: solver
            TYPE (analysisparameters),INTENT(IN)::analysis
            TYPE(macro_GP_DATA),INTENT(INOUT)::GP_DATA
            TYPE(meshparameters),INTENT(IN)::para
            
            IF (analysis%solving_process==0) THEN !full simulation
                CALL solver%PARDISOsolver%finish(para,analysis,GP_DATA)
            ELSE IF (analysis%solving_process==1) THEN !reduced simulation
                CALL solver%LAPACKsolver%finish(analysis)
            END IF
                
        END SUBROUTINE finish_procedure
    
END MODULE class_Solver
