!=================================================================================================================================== ! Copyright (c) 2025 GVEC Contributors, Max Planck Institute for Plasma Physics ! License: MIT !=================================================================================================================================== #include "defines.FPP" !=================================================================================================================================== !> !!# Module **Globals** !! !! Here globally used variables /functions are defined !! !=================================================================================================================================== MODULE MODgvec_Globals #ifndef NOISOENV USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : INPUT_UNIT, OUTPUT_UNIT, ERROR_UNIT #endif USE_MPI ! USE MODgvec_py_abort, ONLY: PY_ABORT IMPLICIT NONE PUBLIC !----------------------------------------------------------------------------------------------------------------------------------- ! Select here the working precision wp !INTEGER, PARAMETER :: wp = selected_real_kind(6,35) !! single precision INTEGER, PARAMETER :: wp = selected_real_kind(15,307) !! double precision !INTEGER, PARAMETER :: wp = selected_real_kind(33,307) !! quadruple precision !----------------------------------------------------------------------------------------------------------------------------------- CHARACTER(LEN=20) :: fmt_sep ='(132("="))' !! formatting of separator line: WRITE(*,fmt_sep) REAL(wp),PARAMETER :: PI =ACOS(-1.0_wp) !! pi parameter REAL(wp),PARAMETER :: TWOPI=2.0_wp*PI !! 2*pi parameter INTEGER :: n_warnings_occured=0 !! for final line in screen output: 0 no warnings occured !----------------------------------------------------------------------------------------------------------------------------------- !for testing LOGICAL :: testdbg=.FALSE. !! for debugging the tests, set true for implementing tests, false to run INTEGER :: testlevel =-1 !! flag for testing routines in code: -1: off INTEGER :: ntestCalled=0 !! counter for called tests INTEGER :: nfailedMsg=0 !! counter for messages on failed tests INTEGER :: testUnit !! unit for out.test file !MPI-------------------------------------------------------------------------------------------------------------------------------- LOGICAL :: MPIRoot=.TRUE. !! flag whether process is MPI root process INTEGER :: myRank=0 !! rank of the MPI task INTEGER :: nRanks=1 !! total number of MPI tasks !----------------------------------------------------------------------------------------------------------------------------------- CHARACTER(LEN=20) :: active_region(5)=(/"", "", "", "", ""/) !! for abort messages, to identify which (sub-)region was currently INTEGER :: iregion=0 !! which active_region to fill INTEGER :: ProgressBar_oldpercent !! for progressBar REAL(wp) :: ProgressBar_starttime !! for progressBar !----------------------------------------------------------------------------------------------------------------------------------- #ifndef NOISOENV INTEGER, PARAMETER :: UNIT_stdIn = input_unit !! Terminal input INTEGER, PARAMETER :: UNIT_stdOut = output_unit !! Terminal output INTEGER, PARAMETER :: UNIT_errOut = error_unit !! For error output #else INTEGER, PARAMETER :: UNIT_stdIn = 5 !! Terminal input INTEGER, PARAMETER :: UNIT_stdOut = 6 !! Terminal output INTEGER, PARAMETER :: UNIT_errOut = 0 !! For error output #endif LOGICAL :: print_backtrace=.TRUE. !! print backtrace on abort if compiled with GNU compiler INTEGER, PARAMETER :: MAXLEN = 4096 !! max length of strings, needed for string handling when compiled with NVHPC INTERFACE reset_subregion MODULE PROCEDURE reset_subregion END INTERFACE INTERFACE enter_subregion MODULE PROCEDURE enter_subregion END INTERFACE INTERFACE exit_subregion MODULE PROCEDURE exit_subregion END INTERFACE INTERFACE Abort MODULE PROCEDURE Abort END INTERFACE ABSTRACT INTERFACE SUBROUTINE RaiseException(ErrorMessage) CHARACTER(LEN=*), INTENT(IN) :: ErrorMessage END SUBROUTINE RaiseException END INTERFACE PROCEDURE(RaiseException), POINTER :: RaiseExceptionPtr => NULL() INTERFACE GetTime MODULE PROCEDURE GetTime END INTERFACE GetTime INTERFACE GetTimeSerial MODULE PROCEDURE GetTimeSerial END INTERFACE GetTimeSerial INTERFACE ProgressBar MODULE PROCEDURE ProgressBar END INTERFACE INTERFACE GETFREEUNIT MODULE PROCEDURE GETFREEUNIT END INTERFACE INTERFACE Eval1DPoly MODULE PROCEDURE Eval1DPoly END INTERFACE INTERFACE CROSS MODULE PROCEDURE CROSS END INTERFACE INTERFACE NORMALIZE MODULE PROCEDURE NORMALIZE END INTERFACE INTERFACE DET33 MODULE PROCEDURE DET33 END INTERFACE INTERFACE INV33 MODULE PROCEDURE INV33 END INTERFACE CONTAINS !================================================================================================================================== !> reset global variables for the subregion output to default !================================================================================================================================== SUBROUTINE reset_subregion ! MODULES IMPLICIT NONE !================================================================================================================================== iregion=0 active_region="" END SUBROUTINE reset_subregion !================================================================================================================================== !> add the current subregion to the active_regions (maximum depth is 5) !! This information is collected uniquely for the abort error message !! !================================================================================================================================== SUBROUTINE enter_subregion(subregion_name) ! MODULES IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES CHARACTER(LEN=*), INTENT(IN) :: subregion_name !---------------------------------------------------------------------------------------------------------------------------------- #if DEBUG CHARACTER(LEN=MAXLEN) :: regions INTEGER :: i #endif !================================================================================================================================== IF(MPIroot)THEN IF(iregion>4) CALL Abort(__STAMP__,& "active subregion reached maximum depth of 5") iregion=iregion+1 active_region(iregion)=subregion_name #if DEBUG regions=active_region(1) DO i=2,iregion regions=TRIM(regions)//"."//TRIM(active_region(i)) END DO SWRITE(Unit_stdOut,'(A)') '==> entering '//TRIM(regions) #endif END IF END SUBROUTINE enter_subregion !================================================================================================================================== !> remove the current subregion from the active subregions !! !================================================================================================================================== SUBROUTINE exit_subregion(subregion_name) ! MODULES IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES CHARACTER(LEN=*), INTENT(IN) :: subregion_name !---------------------------------------------------------------------------------------------------------------------------------- #if DEBUG CHARACTER(LEN=MAXLEN) :: regions INTEGER :: i #endif !================================================================================================================================== IF(MPIroot)THEN #if DEBUG regions=active_region(1) DO i=2,iregion regions=TRIM(regions)//"."//TRIM(active_region(i)) END DO SWRITE(Unit_stdOut,'(A)') '<== exiting '//TRIM(regions) #endif IF(TRIM(subregion_name).NE.TRIM(active_region(iregion))) & CALL Abort(__STAMP__,& "trying to exit subregion '"//TRIM(subregion_name)// & "', but currently active subregion is '"//TRIM(active_region(iregion))//"'") active_region(iregion)="" iregion=iregion-1 END IF END SUBROUTINE exit_subregion !================================================================================================================================== !> Terminate program correctly if an error has occurred (important in MPI mode!). !! Uses a MPI_ABORT which terminates FLUXO if a single proc calls this routine. !! !================================================================================================================================== SUBROUTINE Abort(SourceFile,SourceLine,CompDate,CompTime,ErrorMessage,IntInfo,RealInfo,ErrorCode,TypeInfo) ! MODULES IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES CHARACTER(LEN=*) :: SourceFile !! Source file where error has occurred INTEGER :: SourceLine !! Line in source file CHARACTER(LEN=*) :: CompDate !! Compilation date CHARACTER(LEN=*) :: CompTime !! Compilation time CHARACTER(LEN=*) :: ErrorMessage !! Error message INTEGER,OPTIONAL :: IntInfo !! additional integer value for error message REAL(wp),OPTIONAL :: RealInfo !! additional real value for error message INTEGER,OPTIONAL :: ErrorCode !! used for MPI CHARACTER(LEN=*),OPTIONAL :: TypeInfo !! Error type, default is "RuntimeError". Or e.g. !! "MissingParameterError","InvalidParameterError","FileNotFoundError","InitializationError" !---------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES CHARACTER(LEN=50) :: IntString,RealString,errtype #if MPI INTEGER :: errOut ! Output of MPI_ABORT INTEGER :: signalout ! Output errorcode #endif CHARACTER(LEN=MAXLEN) :: errmsg INTEGER :: i !================================================================================================================================== IntString = "" RealString = "" errtype="RuntimeError" errmsg="" IF(MPIroot)THEN errmsg=TRIM(active_region(1)) DO i=2,iregion errmsg=TRIM(errmsg)//"."//TRIM(active_region(i)) END DO CALL reset_subregion() END IF IF(PRESENT(TypeInfo)) errtype = TRIM(TypeInfo) errmsg=TRIM(errmsg) // " | "//TRIM(errtype) IF (PRESENT(IntInfo)) THEN WRITE(IntString,"(I8)") IntInfo IntString=",IntInfo="//TRIM(IntString) END IF IF (PRESENT(RealInfo)) THEN WRITE(RealString,"(F24.19)") RealInfo RealString=",RealInfo="//TRIM(RealString) END IF errmsg=TRIM(errmsg)//" | "//TRIM(ErrorMessage)//TRIM(IntString)//TRIM(RealString) WRITE(UNIT_stdOut,*) '_____________________________________________________________________________\n', & 'Program abort caused on Proc ',myRank, '\n', & ' in File : ',TRIM(SourceFile),' Line ',SourceLine, '\n', & ' This file was compiled at ',TRIM(CompDate),' ',TRIM(CompTime), '\n', & 'Message: ',TRIM(errmsg) CALL FLUSH(UNIT_stdOut) #if MPI signalout=2 ! MPI_ABORT requires an output error-code /=0 IF(PRESENT(ErrorCode)) signalout=ErrorCode CALL MPI_ABORT(MPI_COMM_WORLD,signalout,errOut) #endif #if GNU IF(print_backtrace) CALL BACKTRACE #endif IF (ASSOCIATED(RaiseExceptionPtr)) THEN CALL RaiseExceptionPtr(errmsg) END IF ERROR STOP 2 END SUBROUTINE Abort !================================================================================================================================== !> Calculates current time (serial / OpenMP /MPI) !! !================================================================================================================================== FUNCTION GetTime() RESULT(t) ! MODULES !$ USE omp_lib IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES REAL(wp) :: t !< output time !---------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES #if MPI LOGICAL :: barr INTEGER :: ierr !================================================================================================================================== CALL MPI_BARRIER(MPI_COMM_WORLD, ierr) ! not possible to 'CALL parBarrier()' because MODgvec_MPI uses MODgvec_Globals! t = MPI_WTIME() #else CALL CPU_TIME(t) !$ t=OMP_GET_WTIME() #endif END FUNCTION GetTime !================================================================================================================================== !> Calculates current time locally on a MPIrank (no MPI Barrier) !! !================================================================================================================================== FUNCTION GetTimeSerial() RESULT(t) ! MODULES !$ USE omp_lib IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES REAL(wp) :: t !< output time !---------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES !================================================================================================================================== CALL CPU_TIME(t) !$ t=OMP_GET_WTIME() END FUNCTION GetTimeSerial !================================================================================================================================== !> Print a progress bar to screen, call either with init=T or init=F !! !================================================================================================================================== SUBROUTINE ProgressBar(iter,n_iter) ! MODULES !$ USE omp_lib IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES INTEGER,INTENT(IN) :: iter,n_iter !! iter ranges from 1...n_iter !---------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES CHARACTER(LEN=8) :: fmtstr INTEGER :: newpercent REAL(wp) :: endTime !================================================================================================================================== IF(.NOT.MPIroot)RETURN IF(iter.LE.0)THEN !INIT ProgressBar_oldpercent=0 ProgressBar_StartTime=GetTimeSerial() WRITE(UNIT_StdOut,'(4X,A,I8)') & '| 10% 20% 30% 40% 50% 60% 70% 80% 90% 100%| ... of ',n_iter WRITE(UNIT_StdOut,'(4X,A1)',ADVANCE='NO')'|' CALL FLUSH(UNIT_stdOut) ELSE newpercent=FLOOR(REAL(iter,wp)/REAL(n_iter,wp)*(100.0_wp+1.0e-12_wp)) WRITE(fmtstr,'(I4)')newpercent-ProgressBar_oldpercent IF(newpercent-ProgressBar_oldpercent.GT.0)THEN WRITE(UNIT_StdOut,'('//TRIM(fmtstr)//'("."))',ADVANCE='NO') CALL FLUSH(UNIT_stdOut) END IF ProgressBar_oldPercent=newPercent IF(newpercent.EQ.100)THEN EndTime=GetTimeSerial() WRITE(Unit_stdOut,'(A3,F8.2,A)') '| [',EndTime-ProgressBar_StartTime,' sec ]' END IF END IF END SUBROUTINE ProgressBar !================================================================================================================================== !> Get unused file unit number !================================================================================================================================== FUNCTION GETFREEUNIT() ! MODULES IMPLICIT NONE !---------------------------------------------------------------------------------------------------------------------------------- ! INPUT/OUTPUT VARIABLES INTEGER :: GetFreeUnit !! File unit number !---------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES LOGICAL :: connected !================================================================================================================================== GetFreeUnit=55 INQUIRE(UNIT=GetFreeUnit, OPENED=connected) IF(connected)THEN DO GetFreeUnit=GetFreeUnit+1 INQUIRE(UNIT=GetFreeUnit, OPENED=connected) IF(.NOT.connected)EXIT END DO END IF END FUNCTION GETFREEUNIT !=================================================================================================================================== !> evalute monomial polynomial c_1+c_2*x+c_3*x^2 ... !! !=================================================================================================================================== PURE FUNCTION Eval1DPoly(nCoefs,Coefs,x) ! MODULES IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES INTEGER, INTENT(IN) :: nCoefs !! number of coefficients REAL(wp), INTENT(IN) :: Coefs(nCoefs) !! coefficients REAL(wp), INTENT(IN) :: x !! evaluation position !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) :: Eval1DPoly !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES INTEGER :: i !=================================================================================================================================== Eval1DPoly=0. DO i=nCoefs,1,-1 Eval1DPoly=Eval1DPoly*x+Coefs(i) END DO END FUNCTION Eval1DPoly !=================================================================================================================================== !> evalute first derivative monomial polynomial (c_1+c_2*x+c_3*x^2) -> (c_2+2*c_3*x+3*c_4*x^2 ... !! !=================================================================================================================================== PURE FUNCTION Eval1DPoly_deriv(nCoefs,Coefs,x) ! MODULES IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES INTEGER, INTENT(IN) :: nCoefs !! number of coefficients REAL(wp), INTENT(IN) :: Coefs(nCoefs) !! coefficients REAL(wp), INTENT(IN) :: x !! evaluation position !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) :: Eval1DPoly_deriv !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES INTEGER :: i !=================================================================================================================================== Eval1DPoly_deriv=0. DO i=nCoefs,2,-1 Eval1DPoly_deriv=Eval1DPoly_deriv*x+REAL(i-1,wp)*Coefs(i) END DO END FUNCTION Eval1DPoly_deriv !=================================================================================================================================== !> normalizes a nDim vector with respect to the eucledian norm !! !=================================================================================================================================== PURE FUNCTION NORMALIZE(v1,nVal) ! MODULES IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES INTEGER,INTENT(IN) :: nVal !! vector size REAL(wp),INTENT(IN) :: v1(nVal) !! vector !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) :: normalize(nVal) !! result, normalized vector !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES !=================================================================================================================================== normalize=v1/SQRT(SUM(v1*v1)) END FUNCTION NORMALIZE !=================================================================================================================================== !> computes the cross product of to 3 dimensional vectors: cross=v1 x v2 !! !=================================================================================================================================== PURE FUNCTION CROSS(v1,v2) ! MODULES IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES REAL(wp),INTENT(IN) :: v1(3) !! first input vector REAL(wp),INTENT(IN) :: v2(3) !! second input vector !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) :: cross(3) !! result v1 x v2 !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES !=================================================================================================================================== cross=(/v1(2)*v2(3)-v1(3)*v2(2),v1(3)*v2(1)-v1(1)*v2(3),v1(1)*v2(2)-v1(2)*v2(1)/) END FUNCTION CROSS !=================================================================================================================================== !> compute determinant of 3x3 matrix !! !=================================================================================================================================== PURE FUNCTION DET33(Mat) ! MODULES IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES REAL(wp),INTENT(IN) :: Mat(3,3) !! input matrix !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) :: DET33 !! determinant of the input matrix !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES !=================================================================================================================================== DET33= ( Mat(1,1) * Mat(2,2) - Mat(1,2) * Mat(2,1) ) * Mat(3,3) & + ( Mat(1,2) * Mat(2,3) - Mat(1,3) * Mat(2,2) ) * Mat(3,1) & + ( Mat(1,3) * Mat(2,1) - Mat(1,1) * Mat(2,3) ) * Mat(3,2) END FUNCTION DET33 !=================================================================================================================================== !> compute inverse of 3x3 matrix, needs sDet=1/det(Mat) !! !=================================================================================================================================== PURE FUNCTION INV33(Mat, Det_in) ! MODULES IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES REAL(wp),INTENT(IN) :: Mat(3,3) !! input matrix REAL(wp),INTENT(IN),OPTIONAL :: Det_in !! determinant of input matrix (otherwise computed here) !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) :: INV33(3,3) !! inverse matrix !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES REAL(wp) :: sDet !=================================================================================================================================== IF(PRESENT(Det_in))THEN sDet=1.0_wp/Det_in ELSE sDet=1.0_wp/DET33(Mat) END IF INV33(1,1) = ( Mat(2,2) * Mat(3,3) - Mat(2,3) * Mat(3,2) ) * sDet INV33(1,2) = ( Mat(1,3) * Mat(3,2) - Mat(1,2) * Mat(3,3) ) * sDet INV33(1,3) = ( Mat(1,2) * Mat(2,3) - Mat(1,3) * Mat(2,2) ) * sDet INV33(2,1) = ( Mat(2,3) * Mat(3,1) - Mat(2,1) * Mat(3,3) ) * sDet INV33(2,2) = ( Mat(1,1) * Mat(3,3) - Mat(1,3) * Mat(3,1) ) * sDet INV33(2,3) = ( Mat(1,3) * Mat(2,1) - Mat(1,1) * Mat(2,3) ) * sDet INV33(3,1) = ( Mat(2,1) * Mat(3,2) - Mat(2,2) * Mat(3,1) ) * sDet INV33(3,2) = ( Mat(1,2) * Mat(3,1) - Mat(1,1) * Mat(3,2) ) * sDet INV33(3,3) = ( Mat(1,1) * Mat(2,2) - Mat(1,2) * Mat(2,1) ) * sDet END FUNCTION INV33 END MODULE MODgvec_Globals