2012-08-09 16:31:53 +05:30
!--------------------------------------------------------------------------------------------------
! $Id$
!--------------------------------------------------------------------------------------------------
2013-09-14 16:29:35 +05:30
!> @author Franz Roters, Max-Planck-Institut für Eisenforschung GmbH
!> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH
!> @author Christoph Kords, Max-Planck-Institut für Eisenforschung GmbH
!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH
2014-05-08 20:25:19 +05:30
!> @author Luv Sharma, Max-Planck-Institut für Eisenforschung GmbH
2013-09-14 16:29:35 +05:30
!> @brief setting precision for real and int type depending on makros "FLOAT" and "INT"
!> @details setting precision for real and int type and for DAMASK_NaN. Definition is made
!! depending on makros "FLOAT" and "INT" defined during compilation
2012-08-09 16:31:53 +05:30
!--------------------------------------------------------------------------------------------------
2012-04-28 16:16:41 +05:30
2012-03-06 20:22:48 +05:30
module prec
2012-08-09 16:31:53 +05:30
2012-03-06 20:22:48 +05:30
implicit none
private
2012-08-28 21:38:17 +05:30
#if (FLOAT==4)
2014-06-27 01:31:38 +05:30
#if defined(Spectral) || defined(FEM)
SPECTRAL SOLVER AND OWN FEM DO NOT SUPPORT SINGLE PRECISION , STOPPING COMPILATION
2012-08-31 01:56:28 +05:30
#endif
2012-08-28 21:38:17 +05:30
integer , parameter , public :: pReal = 4 !< floating point single precition (was selected_real_kind(6,37), number with 6 significant digits, up to 1e+-37)
2012-10-02 22:23:03 +05:30
#ifdef __INTEL_COMPILER
2012-08-28 21:38:17 +05:30
real ( pReal ) , parameter , public :: DAMASK_NaN = Z '7F800001' !< quiet NaN for single precision (from http://www.hpc.unimelb.edu.au/doc/f90lrm/dfum_035.html, copy can be found in documentation/Code/Fortran)
2012-10-02 22:23:03 +05:30
#endif
#ifdef __GFORTRAN__
2012-08-28 21:38:17 +05:30
real ( pReal ) , parameter , public :: DAMASK_NaN = real ( Z '7F800001' , pReal ) !< quiet NaN for single precision (from http://www.hpc.unimelb.edu.au/doc/f90lrm/dfum_035.html, copy can be found in documentation/Code/Fortran)
#endif
#elif (FLOAT==8)
2012-08-29 00:40:54 +05:30
integer , parameter , public :: pReal = 8 !< floating point double precision (was selected_real_kind(15,300), number with 15 significant digits, up to 1e+-300)
2012-10-02 22:23:03 +05:30
#ifdef __INTEL_COMPILER
2012-08-28 21:38:17 +05:30
real ( pReal ) , parameter , public :: DAMASK_NaN = Z '7FF8000000000000' !< quiet NaN for double precision (from http://www.hpc.unimelb.edu.au/doc/f90lrm/dfum_035.html, copy can be found in documentation/Code/Fortran)
2012-10-02 22:23:03 +05:30
#endif
#ifdef __GFORTRAN__
2012-08-28 21:38:17 +05:30
real ( pReal ) , parameter , public :: DAMASK_NaN = real ( Z '7FF8000000000000' , pReal ) !< quiet NaN for double precision (from http://www.hpc.unimelb.edu.au/doc/f90lrm/dfum_035.html, copy can be found in documentation/Code/Fortran)
#endif
2012-08-31 01:56:28 +05:30
#else
2013-09-18 19:29:42 +05:30
NO SUITABLE PRECISION FOR REAL SELECTED , STOPPING COMPILATION
2012-02-13 19:38:07 +05:30
#endif
2007-03-20 19:25:22 +05:30
2012-08-28 21:38:17 +05:30
#if (INT==4)
integer , parameter , public :: pInt = 4 !< integer representation 32 bit (was selected_int_kind(9), number with at least up to +- 1e9)
#elif (INT==8)
integer , parameter , public :: pInt = 8 !< integer representation 64 bit (was selected_int_kind(12), number with at least up to +- 1e12)
2012-08-31 01:56:28 +05:30
#else
2013-09-18 19:29:42 +05:30
NO SUITABLE PRECISION FOR INTEGER SELECTED , STOPPING COMPILATION
2012-08-28 21:38:17 +05:30
#endif
integer , parameter , public :: pLongInt = 8 !< integer representation 64 bit (was selected_int_kind(12), number with at least up to +- 1e12)
2013-09-14 16:29:35 +05:30
real ( pReal ) , parameter , public :: tol_math_check = 1.0e-8_pReal !< tolerance for internal math self-checks (rotation)
2012-08-28 21:38:17 +05:30
2014-10-14 06:03:38 +05:30
type , public :: p_vec !< variable length datatype used for storage of state
real ( pReal ) , dimension ( : ) , allocatable :: p
end type p_vec
type , public :: p_intvec
2014-06-25 21:38:49 +05:30
integer ( pInt ) , dimension ( : ) , allocatable :: p
2013-02-25 18:43:52 +05:30
end type p_intvec
2014-04-15 14:50:38 +05:30
!http://stackoverflow.com/questions/3948210/can-i-have-a-pointer-to-an-item-in-an-allocatable-array
type , public :: tState
2014-12-09 17:42:53 +05:30
integer ( pInt ) :: &
sizeState = 0_pInt , & !< size of state
sizeDotState = 0_pInt , & !< size of dot state, i.e. parts of the state that are integrated
sizePostResults = 0_pInt !< size of output data
logical :: &
nonlocal = . false . !< absolute tolerance for state integration
real ( pReal ) , allocatable , dimension ( : ) :: &
atolState
real ( pReal ) , pointer , dimension ( : , : ) :: &
state , & !< state
dotState !< state rate
real ( pReal ) , allocatable , dimension ( : , : ) :: &
state0 , &
partionedState0 , &
subState0 , &
state_backup , &
deltaState , &
previousDotState , & !< state rate of previous xxxx
previousDotState2 , & !< state rate two xxxx ago
dotState_backup , & !< backup of state rate
RK4dotState
real ( pReal ) , allocatable , dimension ( : , : , : ) :: &
RKCK45dotState
end type
type , extends ( tState ) , public :: tPlasticState
integer ( pInt ) :: &
nSlip = 0_pInt , &
nTwin = 0_pInt , &
nTrans = 0_pInt
real ( pReal ) , pointer , dimension ( : , : ) :: &
slipRate , & !< slip rate
accumulatedSlip !< accumulated plastic slip
2014-04-15 14:50:38 +05:30
end type
2014-09-19 23:29:06 +05:30
type , public :: tFieldData
2014-12-09 17:42:53 +05:30
integer ( pInt ) :: &
sizeField = 0_pInt , &
sizePostResults = 0_pInt
real ( pReal ) , allocatable , dimension ( : , : ) :: &
field !< field data
2014-08-21 23:18:20 +05:30
end type
2014-04-15 14:50:38 +05:30
2013-02-11 15:14:17 +05:30
public :: &
prec_init
2012-03-06 20:22:48 +05:30
contains
2012-08-28 21:38:17 +05:30
2012-08-09 16:31:53 +05:30
!--------------------------------------------------------------------------------------------------
!> @brief reporting precision and checking if DAMASK_NaN is set correctly
!--------------------------------------------------------------------------------------------------
2012-03-06 20:22:48 +05:30
subroutine prec_init
2012-08-09 16:31:53 +05:30
use , intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment)
2012-03-06 20:22:48 +05:30
implicit none
2014-10-10 18:38:34 +05:30
integer ( pInt ) :: worldrank = 0_pInt
#ifdef PETSc
2015-03-18 22:48:43 +05:30
#include <petsc-finclude/petscsys.h>
2014-10-10 01:53:06 +05:30
PetscErrorCode :: ierr
#endif
2013-02-13 23:24:56 +05:30
external :: &
2015-03-29 18:24:13 +05:30
quit , &
MPI_Comm_rank , &
MPI_Abort
2013-02-13 23:24:56 +05:30
2014-10-10 18:38:34 +05:30
#ifdef PETSc
2014-10-10 01:53:06 +05:30
call MPI_Comm_rank ( PETSC_COMM_WORLD , worldrank , ierr ) ; CHKERRQ ( ierr )
2014-10-10 18:38:34 +05:30
#endif
mainProcess : if ( worldrank == 0 ) then
write ( 6 , '(/,a)' ) ' <<<+- prec init -+>>>'
write ( 6 , '(a)' ) ' $Id$'
2012-02-01 00:48:55 +05:30
#include "compilation_info.f90"
2014-10-10 18:38:34 +05:30
write ( 6 , '(a,i3)' ) ' Bytes for pReal: ' , pReal
write ( 6 , '(a,i3)' ) ' Bytes for pInt: ' , pInt
write ( 6 , '(a,i3)' ) ' Bytes for pLongInt: ' , pLongInt
write ( 6 , '(a,e10.3)' ) ' NaN: ' , DAMASK_NaN
write ( 6 , '(a,l3,/)' ) ' NaN /= NaN: ' , DAMASK_NaN / = DAMASK_NaN
endif mainProcess
2014-04-15 14:50:38 +05:30
2012-10-04 19:52:39 +05:30
if ( DAMASK_NaN == DAMASK_NaN ) call quit ( 9000 )
2009-08-31 20:39:15 +05:30
2012-03-06 20:22:48 +05:30
end subroutine prec_init
2009-08-31 20:39:15 +05:30
2012-03-06 20:22:48 +05:30
end module prec