DAMASK (except nonlocal) now sees and interacts with FEM solver
This commit is contained in:
Pratheek Shanthraj
parent
9a285a1e9f
commit
639ca89133
|
|
@ -71,6 +71,10 @@ subroutine CPFEM_initAll(temperature,el,ip)
|
|||
use IO, only: &
|
||||
IO_init
|
||||
use DAMASK_interface
|
||||
#ifdef FEM
|
||||
use FEZoo, only: &
|
||||
FEZoo_init
|
||||
#endif
|
||||
|
||||
implicit none
|
||||
integer(pInt), intent(in) :: el, & ! FE el number
|
||||
|
|
@ -81,9 +85,15 @@ subroutine CPFEM_initAll(temperature,el,ip)
|
|||
if (.not. CPFEM_init_done) then
|
||||
#ifdef Spectral
|
||||
call DAMASK_interface_init() ! Spectral solver is interfacing to commandline
|
||||
#endif
|
||||
#ifdef FEM
|
||||
call DAMASK_interface_init() ! Spectral solver is interfacing to commandline
|
||||
#endif
|
||||
call prec_init
|
||||
call IO_init
|
||||
#ifdef FEM
|
||||
call FEZoo_init() ! Spectral solver is interfacing to commandline
|
||||
#endif
|
||||
call numerics_init
|
||||
call debug_init
|
||||
call math_init
|
||||
|
|
@ -96,7 +106,9 @@ subroutine CPFEM_initAll(temperature,el,ip)
|
|||
call homogenization_init
|
||||
call CPFEM_init
|
||||
#ifndef Spectral
|
||||
#ifndef FEM
|
||||
call DAMASK_interface_init() ! Spectral solver init is already done
|
||||
#endif
|
||||
#endif
|
||||
CPFEM_init_done = .true.
|
||||
endif
|
||||
|
|
|
|||
|
|
@ -67,8 +67,10 @@ subroutine FE_init
|
|||
IO_intValue, &
|
||||
IO_lc, &
|
||||
#ifndef Spectral
|
||||
#ifndef FEM
|
||||
IO_open_inputFile, &
|
||||
IO_open_logFile, &
|
||||
#endif
|
||||
#endif
|
||||
IO_warning, &
|
||||
IO_timeStamp
|
||||
|
|
@ -76,6 +78,7 @@ subroutine FE_init
|
|||
|
||||
implicit none
|
||||
#ifndef Spectral
|
||||
#ifndef FEM
|
||||
integer(pInt), parameter :: &
|
||||
FILEUNIT = 222_pInt, &
|
||||
maxNchunks = 6_pInt
|
||||
|
|
@ -83,6 +86,7 @@ subroutine FE_init
|
|||
character(len=64) :: tag
|
||||
character(len=1024) :: line
|
||||
integer(pInt), dimension(1_pInt+2_pInt*maxNchunks) :: positions
|
||||
#endif
|
||||
#endif
|
||||
|
||||
write(6,'(/,a)') ' <<<+- FEsolving init -+>>>'
|
||||
|
|
@ -98,6 +102,13 @@ subroutine FE_init
|
|||
restartInc = 1_pInt
|
||||
endif
|
||||
restartRead = restartInc > 1_pInt ! only read in if "true" restart requested
|
||||
#elif defined FEM
|
||||
restartInc = FEMRestartInc
|
||||
if(restartInc <= 0_pInt) then
|
||||
call IO_warning(warning_ID=34_pInt)
|
||||
restartInc = 1_pInt
|
||||
endif
|
||||
restartRead = restartInc > 1_pInt
|
||||
#else
|
||||
call IO_open_inputFile(FILEUNIT,modelName)
|
||||
rewind(FILEUNIT)
|
||||
|
|
|
|||
|
|
@ -73,10 +73,12 @@ module IO
|
|||
IO_intOut, &
|
||||
IO_timeStamp
|
||||
#ifndef Spectral
|
||||
#ifndef FEM
|
||||
public :: &
|
||||
IO_open_inputFile, &
|
||||
IO_open_logFile
|
||||
#endif
|
||||
#endif
|
||||
#ifdef Abaqus
|
||||
public :: &
|
||||
IO_abaqus_hasNoPart
|
||||
|
|
@ -302,6 +304,7 @@ end function IO_open_JobFile_stat
|
|||
|
||||
|
||||
#ifndef Spectral
|
||||
#ifndef FEM
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
!> @brief opens FEM input file for reading located in current working directory to given unit
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
|
|
@ -367,6 +370,7 @@ subroutine IO_open_logFile(fileUnit)
|
|||
|
||||
end subroutine IO_open_logFile
|
||||
#endif
|
||||
#endif
|
||||
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
|
|
|
|||
|
|
@ -29,6 +29,9 @@ module constitutive
|
|||
constitutive_aTolState !< pointer array to absolute state tolerance
|
||||
type(p_vec), public, dimension(:,:,:,:), allocatable :: &
|
||||
constitutive_RKCK45dotState !< pointer array to evolution of microstructure used by Cash-Karp Runge-Kutta method
|
||||
real(pReal), public, dimension(:,:,:), allocatable :: &
|
||||
constitutive_localDamage, &
|
||||
constitutive_gradientDamage
|
||||
integer(pInt), public, dimension(:,:,:), allocatable :: &
|
||||
constitutive_sizeDotState, & !< size of dotState array
|
||||
constitutive_sizeState, & !< size of state array per grain
|
||||
|
|
@ -240,6 +243,8 @@ subroutine constitutive_init
|
|||
#ifndef NEWSTATE
|
||||
! lumped into new state
|
||||
allocate(constitutive_state0(cMax,iMax,eMax))
|
||||
allocate(constitutive_localDamage(cMax,iMax,eMax)); constitutive_localDamage = 1.0_pReal
|
||||
allocate(constitutive_gradientDamage(cMax,iMax,eMax)); constitutive_gradientDamage = 1.0_pReal
|
||||
allocate(constitutive_partionedState0(cMax,iMax,eMax))
|
||||
allocate(constitutive_subState0(cMax,iMax,eMax))
|
||||
allocate(constitutive_state(cMax,iMax,eMax))
|
||||
|
|
@ -775,7 +780,7 @@ use math, only : &
|
|||
C = math_Mandel66to3333(constitutive_homogenizedC(ipc,ip,el))
|
||||
|
||||
FeT = math_transpose33(Fe)
|
||||
T = 0.5_pReal * math_mul3333xx33(C, math_mul33x33(FeT,Fe)-MATH_I3)
|
||||
T = 0.5_pReal*math_mul3333xx33(C,math_mul33x33(FeT,Fe)-MATH_I3)*constitutive_gradientDamage(ipc,ip,el)
|
||||
|
||||
dT_dFe = 0.0_pReal
|
||||
forall (i=1_pInt:3_pInt, j=1_pInt:3_pInt, k=1_pInt:3_pInt, l=1_pInt:3_pInt) &
|
||||
|
|
|
|||
|
|
@ -551,7 +551,9 @@ subroutine crystallite_stressAndItsTangent(updateJaco,rate_sensitivity)
|
|||
#else
|
||||
mappingConstitutive, &
|
||||
#endif
|
||||
constitutive_TandItsTangent
|
||||
constitutive_TandItsTangent, &
|
||||
constitutive_localDamage, &
|
||||
constitutive_gradientDamage
|
||||
|
||||
implicit none
|
||||
logical, intent(in) :: &
|
||||
|
|
@ -1371,6 +1373,14 @@ subroutine crystallite_stressAndItsTangent(updateJaco,rate_sensitivity)
|
|||
crystallite_heat(g,i,e) = 0.98_pReal* &
|
||||
abs(math_mul33xx33(math_Mandel6to33(crystallite_Tstar_v(1:6,g,i,e)), &
|
||||
crystallite_Lp(1:3,1:3,g,i,e)))
|
||||
constitutive_localDamage(g,i,e) = &
|
||||
1.0_pReal* &
|
||||
sum(math_Mandel6to33(crystallite_Tstar_v(1:6,g,i,e)/constitutive_gradientDamage(g,i,e))* &
|
||||
(math_mul33x33(math_transpose33(crystallite_Fe(1:3,1:3,g,i,e)), &
|
||||
crystallite_Fe(1:3,1:3,g,i,e))-math_I3))/4.0_pReal + &
|
||||
0.0_pReal* &
|
||||
sum(abs(math_mul33x33(math_transpose33(crystallite_Fp(1:3,1:3,g,i,e)), &
|
||||
crystallite_Fp(1:3,1:3,g,i,e))-math_I3)/2.0_pReal)
|
||||
enddo
|
||||
enddo
|
||||
enddo elementLooping12
|
||||
|
|
|
|||
|
|
@ -11,6 +11,9 @@ end module libs
|
|||
#ifdef Spectral
|
||||
#include "../lib/kdtree2.f90"
|
||||
#endif
|
||||
#ifdef FEM
|
||||
#include "../lib/kdtree2.f90"
|
||||
#endif
|
||||
#include "../lib/IR_Precision.f90"
|
||||
#include "../lib/Lib_Base64.f90"
|
||||
#include "../lib/Lib_VTK_IO.f90"
|
||||
|
|
|
|||
|
|
@ -12,6 +12,9 @@ module numerics
|
|||
|
||||
implicit none
|
||||
private
|
||||
#ifdef FEM
|
||||
#include <finclude/petsc.h90>
|
||||
#endif
|
||||
character(len=64), parameter, private :: &
|
||||
numerics_CONFIGFILE = 'numerics.config' !< name of configuration file
|
||||
|
||||
|
|
@ -89,6 +92,7 @@ module numerics
|
|||
fftw_planner_flag = 32_pInt, & !< conversion of fftw_plan_mode to integer, basically what is usually done in the include file of fftw
|
||||
itmax = 250_pInt, & !< maximum number of iterations
|
||||
itmin = 2_pInt, & !< minimum number of iterations
|
||||
stagItMax = 25_pInt, & !< maximum number of iterations
|
||||
maxCutBack = 3_pInt, & !< max number of cut backs
|
||||
continueCalculation = 0_pInt, & !< 0: exit if BVP solver does not converge, 1: continue calculation if BVP solver does not converge
|
||||
divergence_correction = 2_pInt !< correct divergence calculation in fourier space 0: no correction, 1: size scaled to 1, 2: size scaled to Npoints
|
||||
|
|
@ -97,6 +101,33 @@ module numerics
|
|||
update_gamma = .false. !< update gamma operator with current stiffness, Default .false.: use initial stiffness
|
||||
#endif
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! FEM parameters:
|
||||
#ifdef FEM
|
||||
real(pReal), protected, public :: &
|
||||
err_struct_tolAbs = 1.0e-10_pReal, & !< absolute tolerance for equilibrium
|
||||
err_struct_tolRel = 5.0e-4_pReal, & !< relative tolerance for equilibrium
|
||||
err_thermal_tol = 1.0e-6_pReal, &
|
||||
err_damage_tol = 1.0e-4_pReal
|
||||
character(len=1024), protected, public :: &
|
||||
petsc_optionsFEM = '-snes_type ngmres &
|
||||
&-snes_ngmres_anderson '
|
||||
integer(pInt), protected, public :: &
|
||||
itmaxFEM = 250_pInt, & !< maximum number of iterations
|
||||
itminFEM = 2_pInt, & !< minimum number of iterations
|
||||
maxCutBackFEM = 3_pInt, & !< max number of cut backs
|
||||
integrationOrder = 1_pInt, &
|
||||
structOrder = 2_pInt, &
|
||||
thermalOrder = 2_pInt, &
|
||||
damageOrder = 2_pInt, &
|
||||
worldrank = 0_pInt, &
|
||||
worldsize = 0_pInt
|
||||
logical, protected, public :: &
|
||||
structThermalCoupling = .false., &
|
||||
structDamageCoupling = .false., &
|
||||
thermalDamageCoupling = .false.
|
||||
#endif
|
||||
|
||||
public :: numerics_init
|
||||
|
||||
contains
|
||||
|
|
@ -124,15 +155,20 @@ subroutine numerics_init
|
|||
|
||||
#ifdef Spectral
|
||||
!$ use OMP_LIB, only: omp_set_num_threads ! Use the standard conforming module file for omp if using the spectral solver
|
||||
#endif
|
||||
#ifdef FEM
|
||||
!$ use OMP_LIB, only: omp_set_num_threads ! Use the standard conforming module file for omp if using the spectral solver
|
||||
#endif
|
||||
implicit none
|
||||
#ifndef Spectral
|
||||
#ifndef FEM
|
||||
!$ include "omp_lib.h" ! use the not F90 standard conforming include file to prevent crashes with some versions of MSC.Marc
|
||||
#endif
|
||||
#endif
|
||||
integer(pInt), parameter :: FILEUNIT = 300_pInt ,&
|
||||
maxNchunks = 2_pInt
|
||||
!$ integer :: gotDAMASK_NUM_THREADS = 1
|
||||
integer :: i ! no pInt
|
||||
integer :: i, ierr ! no pInt
|
||||
integer(pInt), dimension(1+2*maxNchunks) :: positions
|
||||
character(len=65536) :: &
|
||||
tag ,&
|
||||
|
|
@ -154,6 +190,11 @@ subroutine numerics_init
|
|||
!$ endif
|
||||
!$ call omp_set_num_threads(DAMASK_NumThreadsInt) ! set number of threads for parallel execution
|
||||
|
||||
#ifdef FEM
|
||||
call MPI_Comm_rank(PETSC_COMM_WORLD,worldrank,ierr);CHKERRQ(ierr)
|
||||
call MPI_Comm_size(PETSC_COMM_WORLD,worldsize,ierr);CHKERRQ(ierr)
|
||||
#endif
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! try to open the config file
|
||||
fileExists: if(IO_open_file_stat(FILEUNIT,numerics_configFile)) then
|
||||
|
|
@ -274,6 +315,8 @@ subroutine numerics_init
|
|||
itmax = IO_intValue(line,positions,2_pInt)
|
||||
case ('itmin')
|
||||
itmin = IO_intValue(line,positions,2_pInt)
|
||||
case ('stagitmax')
|
||||
stagItMax = IO_intValue(line,positions,2_pInt)
|
||||
case ('maxcutback')
|
||||
maxCutBack = IO_intValue(line,positions,2_pInt)
|
||||
case ('continuecalculation')
|
||||
|
|
@ -312,12 +355,53 @@ subroutine numerics_init
|
|||
#endif
|
||||
#ifndef Spectral
|
||||
case ('err_div_tolabs','err_div_tolrel','err_stress_tolrel','err_stress_tolabs',& ! found spectral parameter for FEM build
|
||||
'itmax', 'itmin','memory_efficient','fftw_timelimit','fftw_plan_mode', &
|
||||
'itmax', 'itmin','stagitmax','memory_efficient','fftw_timelimit','fftw_plan_mode', &
|
||||
'divergence_correction','update_gamma','myfilter', &
|
||||
'err_curl_tolabs','err_curl_tolrel', &
|
||||
'maxcutback','polaralpha','polarbeta')
|
||||
call IO_warning(40_pInt,ext_msg=tag)
|
||||
#endif
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! FEM parameters
|
||||
#ifdef FEM
|
||||
case ('err_struct_tolabs')
|
||||
err_struct_tolAbs = IO_floatValue(line,positions,2_pInt)
|
||||
case ('err_struct_tolrel')
|
||||
err_struct_tolRel = IO_floatValue(line,positions,2_pInt)
|
||||
case ('err_thermal_tol')
|
||||
err_thermal_tol = IO_floatValue(line,positions,2_pInt)
|
||||
case ('err_damage_tol')
|
||||
err_damage_tol = IO_floatValue(line,positions,2_pInt)
|
||||
case ('itmaxfem')
|
||||
itmaxFEM = IO_intValue(line,positions,2_pInt)
|
||||
case ('itminfem')
|
||||
itminFEM = IO_intValue(line,positions,2_pInt)
|
||||
case ('maxcutbackfem')
|
||||
maxCutBackFEM = IO_intValue(line,positions,2_pInt)
|
||||
case ('integrationorder')
|
||||
integrationorder = IO_intValue(line,positions,2_pInt)
|
||||
case ('structorder')
|
||||
structorder = IO_intValue(line,positions,2_pInt)
|
||||
case ('thermalorder')
|
||||
thermalorder = IO_intValue(line,positions,2_pInt)
|
||||
case ('damageorder')
|
||||
damageorder = IO_intValue(line,positions,2_pInt)
|
||||
case ('petsc_optionsfem')
|
||||
petsc_optionsFEM = trim(line(positions(4):))
|
||||
case ('structthermalcoupling')
|
||||
structThermalCoupling = IO_intValue(line,positions,2_pInt) > 0_pInt
|
||||
case ('structdamagecoupling')
|
||||
structDamageCoupling = IO_intValue(line,positions,2_pInt) > 0_pInt
|
||||
case ('thermaldamagecoupling')
|
||||
thermalDamageCoupling = IO_intValue(line,positions,2_pInt) > 0_pInt
|
||||
#endif
|
||||
#ifndef FEM
|
||||
case ('err_struct_tolabs','err_struct_tolrel','err_thermal_tol','err_damage_tol',& ! found spectral parameter for FEM build
|
||||
'itmaxfem', 'itminfem','maxcutbackfem','integrationorder','structorder','thermalorder', &
|
||||
'damageorder','petsc_optionsfem','structthermalcoupling','structdamagecoupling','thermaldamagecoupling')
|
||||
call IO_warning(40_pInt,ext_msg=tag)
|
||||
#endif
|
||||
case default ! found unknown keyword
|
||||
call IO_error(300_pInt,ext_msg=tag)
|
||||
endselect
|
||||
|
|
@ -405,6 +489,7 @@ subroutine numerics_init
|
|||
#ifdef Spectral
|
||||
write(6,'(a24,1x,i8)') ' itmax: ',itmax
|
||||
write(6,'(a24,1x,i8)') ' itmin: ',itmin
|
||||
write(6,'(a24,1x,i8)') ' stagItMax: ',stagItMax
|
||||
write(6,'(a24,1x,i8)') ' maxCutBack: ',maxCutBack
|
||||
write(6,'(a24,1x,i8)') ' continueCalculation: ',continueCalculation
|
||||
write(6,'(a24,1x,L8)') ' memory_efficient: ',memory_efficient
|
||||
|
|
@ -432,6 +517,26 @@ subroutine numerics_init
|
|||
#endif
|
||||
#endif
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! spectral parameters
|
||||
#ifdef FEM
|
||||
write(6,'(a24,1x,i8)') ' itmaxFEM: ',itmaxFEM
|
||||
write(6,'(a24,1x,i8)') ' itminFEM: ',itminFEM
|
||||
write(6,'(a24,1x,i8)') ' maxCutBackFEM: ',maxCutBackFEM
|
||||
write(6,'(a24,1x,i8)') ' integrationOrder: ',integrationOrder
|
||||
write(6,'(a24,1x,i8)') ' structOrder: ',structOrder
|
||||
write(6,'(a24,1x,i8)') ' thermalOrder: ',thermalOrder
|
||||
write(6,'(a24,1x,i8)') ' damageOrder: ',damageOrder
|
||||
write(6,'(a24,1x,es8.1)') ' err_struct_tolAbs: ',err_struct_tolAbs
|
||||
write(6,'(a24,1x,es8.1)') ' err_struct_tolRel: ',err_struct_tolRel
|
||||
write(6,'(a24,1x,es8.1)') ' err_thermal_tol: ',err_thermal_tol
|
||||
write(6,'(a24,1x,es8.1)') ' err_damage_tol: ',err_damage_tol
|
||||
write(6,'(a24,1x,L8)') ' structThermalCoupling: ',structThermalCoupling
|
||||
write(6,'(a24,1x,L8)') ' structDamageCoupling: ',structDamageCoupling
|
||||
write(6,'(a24,1x,L8)') ' thermalDamageCoupling: ',thermalDamageCoupling
|
||||
write(6,'(a24,1x,a)') ' PETSc_optionsFEM: ',trim(petsc_optionsFEM)
|
||||
#endif
|
||||
|
||||
!--------------------------------------------------------------------------------------------------
|
||||
! sanity checks
|
||||
if (relevantStrain <= 0.0_pReal) call IO_error(301_pInt,ext_msg='relevantStrain')
|
||||
|
|
@ -473,7 +578,7 @@ subroutine numerics_init
|
|||
if (volDiscrPow_RGC <= 0.0_pReal) call IO_error(301_pInt,ext_msg='volDiscrPw_RGC')
|
||||
#ifdef Spectral
|
||||
if (itmax <= 1_pInt) call IO_error(301_pInt,ext_msg='itmax')
|
||||
if (itmin > itmax .or. itmin < 1_pInt) call IO_error(301_pInt,ext_msg='itmin')
|
||||
if (itmin > itmax .or. itmin < 0_pInt) call IO_error(301_pInt,ext_msg='itmin')
|
||||
if (continueCalculation /= 0_pInt .and. &
|
||||
continueCalculation /= 1_pInt) call IO_error(301_pInt,ext_msg='continueCalculation')
|
||||
if (divergence_correction < 0_pInt .or. &
|
||||
|
|
@ -494,6 +599,16 @@ subroutine numerics_init
|
|||
polarBeta > 2.0_pReal) call IO_error(301_pInt,ext_msg='polarBeta')
|
||||
#endif
|
||||
#endif
|
||||
#ifdef FEM
|
||||
if (itmaxFEM <= 1_pInt) call IO_error(301_pInt,ext_msg='itmaxFEM')
|
||||
if (itminFEM > itmaxFEM .or. &
|
||||
itminFEM < 0_pInt) call IO_error(301_pInt,ext_msg='itminFEM')
|
||||
if (maxCutBackFEM < 0_pInt) call IO_error(301_pInt,ext_msg='maxCutBackFEM')
|
||||
if (err_struct_tolRel <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_struct_tolRel')
|
||||
if (err_struct_tolAbs <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_struct_tolAbs')
|
||||
if (err_thermal_tol <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_thermal_tol')
|
||||
if (err_damage_tol <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_damage_tol')
|
||||
#endif
|
||||
|
||||
end subroutine numerics_init
|
||||
|
||||
|
|
|
|||
|
|
@ -18,6 +18,9 @@ module prec
|
|||
#if (FLOAT==4)
|
||||
#ifdef Spectral
|
||||
SPECTRAL SOLVER DOES NOT SUPPORT SINGLE PRECISION, STOPPING COMPILATION
|
||||
#endif
|
||||
#ifdef FEM
|
||||
SPECTRAL SOLVER DOES NOT SUPPORT SINGLE PRECISION, STOPPING COMPILATION
|
||||
#endif
|
||||
integer, parameter, public :: pReal = 4 !< floating point single precition (was selected_real_kind(6,37), number with 6 significant digits, up to 1e+-37)
|
||||
#ifdef __INTEL_COMPILER
|
||||
|
|
|
|||
Loading…
Reference in New Issue