223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

unified naming scheme, fixed a bug in the new basic scheme, and added a statistic file similar to abaqus to store the information on needed cut backs and iterations for each inc

This commit is contained in:
Martin Diehl 2012-10-19 08:44:21 +00:00
parent 6230dacbac
commit b9f97ba5da
6 changed files with 115 additions and 82 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

84
code/DAMASK_spectral_Driver.f90 → code/DAMASK_spectral_driver.f90
View File

@ -71,22 +71,22 @@ program DAMASK_spectral_Driver
implicit none implicit none
type tLoadCase type tLoadCase
real(pReal), dimension (3,3) :: rotation = math_I3 ! rotation of BC real(pReal), dimension (3,3) :: rotation = math_I3 !< rotation of BC
type(tBoundaryCondition) :: P, & ! stress BC type(tBoundaryCondition) :: P, & !< stress BC
deformation ! deformation BC (Fdot or L) deformation !< deformation BC (Fdot or L)
real(pReal) :: time = 0.0_pReal, & ! length of increment real(pReal) :: time = 0.0_pReal, & !< length of increment
temperature = 300.0_pReal ! isothermal starting conditions temperature = 300.0_pReal !< isothermal starting conditions
integer(pInt) :: incs = 0_pInt, & ! number of increments integer(pInt) :: incs = 0_pInt, & !< number of increments
outputfrequency = 1_pInt, & ! frequency of result writes outputfrequency = 1_pInt, & !< frequency of result writes
restartfrequency = 0_pInt, & ! frequency of restart writes restartfrequency = 0_pInt, & !< frequency of restart writes
logscale = 0_pInt ! linear/logaritmic time inc flag logscale = 0_pInt !< linear/logaritmic time inc flag
logical :: followFormerTrajectory = .true. ! follow trajectory of former loadcase logical :: followFormerTrajectory = .true. !< follow trajectory of former loadcase
end type tLoadCase end type tLoadCase
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! variables related to information from load case and geom file ! variables related to information from load case and geom file
real(pReal), dimension(9) :: temp_valueVector !> temporarily from loadcase file when reading in tensors real(pReal), dimension(9) :: temp_valueVector !< temporarily from loadcase file when reading in tensors
logical, dimension(9) :: temp_maskVector !> temporarily from loadcase file when reading in tensors logical, dimension(9) :: temp_maskVector !< temporarily from loadcase file when reading in tensors
integer(pInt), parameter :: maxNchunksLoadcase = (1_pInt + 9_pInt)*3_pInt +& ! deformation, rotation, and stress integer(pInt), parameter :: maxNchunksLoadcase = (1_pInt + 9_pInt)*3_pInt +& ! deformation, rotation, and stress
(1_pInt + 1_pInt)*5_pInt +& ! time, (log)incs, temp, restartfrequency, and outputfrequency (1_pInt + 1_pInt)*5_pInt +& ! time, (log)incs, temp, restartfrequency, and outputfrequency
1_pInt ! dropguessing 1_pInt ! dropguessing
@ -104,13 +104,14 @@ program DAMASK_spectral_Driver
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! loop variables, convergence etc. ! loop variables, convergence etc.
real(pReal), dimension(3,3), parameter :: ones = 1.0_pReal, zeroes = 0.0_pReal real(pReal), dimension(3,3), parameter :: ones = 1.0_pReal, zeros = 0.0_pReal
real(pReal) :: time = 0.0_pReal, time0 = 0.0_pReal, timeinc = 1.0_pReal, timeinc_old = 0.0_pReal ! elapsed time, begin of interval, time interval, previous time interval real(pReal) :: time = 0.0_pReal, time0 = 0.0_pReal, timeinc = 1.0_pReal, timeinc_old = 0.0_pReal ! elapsed time, begin of interval, time interval, previous time interval
real(pReal) :: guessmode real(pReal) :: guessmode
integer(pInt) :: i, j, k, l, errorID, cutBackLevel = 0_pInt, stepFraction = 0_pInt integer(pInt) :: i, j, k, l, errorID, cutBackLevel = 0_pInt, stepFraction = 0_pInt
integer(pInt) :: currentLoadcase = 0_pInt, inc, & integer(pInt) :: currentLoadcase = 0_pInt, inc, &
totalIncsCounter = 0_pInt,& totalIncsCounter = 0_pInt,&
notConvergedCounter = 0_pInt, convergedCounter = 0_pInt notConvergedCounter = 0_pInt, convergedCounter = 0_pInt, &
resUnit = 0_pInt, statUnit = 0_pInt
character(len=6) :: loadcase_string character(len=6) :: loadcase_string
character(len=1024) :: incInfo character(len=1024) :: incInfo
type(tLoadCase), allocatable, dimension(:) :: loadCases type(tLoadCase), allocatable, dimension(:) :: loadCases
@ -121,7 +122,7 @@ program DAMASK_spectral_Driver
call CPFEM_initAll(temperature = 300.0_pReal, element = 1_pInt, IP= 1_pInt) call CPFEM_initAll(temperature = 300.0_pReal, element = 1_pInt, IP= 1_pInt)
write(6,'(a)') '' write(6,'(a)') ''
write(6,'(a)') ' <<<+- DAMASK_spectral_Driver init -+>>>' write(6,'(a)') ' <<<+- DAMASK_spectral_driver init -+>>>'
write(6,'(a)') ' $Id$' write(6,'(a)') ' $Id$'
#include "compilation_info.f90" #include "compilation_info.f90"
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -175,7 +176,7 @@ program DAMASK_spectral_Driver
if (temp_maskVector(j)) temp_valueVector(j) = IO_floatValue(line,positions,i+j) if (temp_maskVector(j)) temp_valueVector(j) = IO_floatValue(line,positions,i+j)
enddo enddo
loadCases(currentLoadCase)%deformation%maskLogical = transpose(reshape(temp_maskVector,[ 3,3])) loadCases(currentLoadCase)%deformation%maskLogical = transpose(reshape(temp_maskVector,[ 3,3]))
loadCases(currentLoadCase)%deformation%maskFloat = merge(ones,zeroes,& loadCases(currentLoadCase)%deformation%maskFloat = merge(ones,zeros,&
loadCases(currentLoadCase)%deformation%maskLogical) loadCases(currentLoadCase)%deformation%maskLogical)
loadCases(currentLoadCase)%deformation%values = math_plain9to33(temp_valueVector) loadCases(currentLoadCase)%deformation%values = math_plain9to33(temp_valueVector)
case('p','pk1','piolakirchhoff','stress') case('p','pk1','piolakirchhoff','stress')
@ -185,7 +186,7 @@ program DAMASK_spectral_Driver
if (temp_maskVector(j)) temp_valueVector(j) = IO_floatValue(line,positions,i+j) if (temp_maskVector(j)) temp_valueVector(j) = IO_floatValue(line,positions,i+j)
enddo enddo
loadCases(currentLoadCase)%P%maskLogical = transpose(reshape(temp_maskVector,[ 3,3])) loadCases(currentLoadCase)%P%maskLogical = transpose(reshape(temp_maskVector,[ 3,3]))
loadCases(currentLoadCase)%P%maskFloat = merge(ones,zeroes,& loadCases(currentLoadCase)%P%maskFloat = merge(ones,zeros,&
loadCases(currentLoadCase)%P%maskLogical) loadCases(currentLoadCase)%P%maskLogical)
loadCases(currentLoadCase)%P%values = math_plain9to33(temp_valueVector) loadCases(currentLoadCase)%P%values = math_plain9to33(temp_valueVector)
case('t','time','delta') ! increment time case('t','time','delta') ! increment time
@ -256,7 +257,6 @@ program DAMASK_spectral_Driver
write(6,'(2x,a,i5)') 'increments: ', loadCases(currentLoadCase)%incs write(6,'(2x,a,i5)') 'increments: ', loadCases(currentLoadCase)%incs
write(6,'(2x,a,i5)') 'output frequency: ', loadCases(currentLoadCase)%outputfrequency write(6,'(2x,a,i5)') 'output frequency: ', loadCases(currentLoadCase)%outputfrequency
write(6,'(2x,a,i5,/)') 'restart frequency: ', loadCases(currentLoadCase)%restartfrequency write(6,'(2x,a,i5,/)') 'restart frequency: ', loadCases(currentLoadCase)%restartfrequency
if (any(loadCases(currentLoadCase)%P%maskLogical .eqv. loadCases(currentLoadCase)%deformation%maskLogical)) errorID = 831_pInt ! exclusive or masking only if (any(loadCases(currentLoadCase)%P%maskLogical .eqv. loadCases(currentLoadCase)%deformation%maskLogical)) errorID = 831_pInt ! exclusive or masking only
if (any(loadCases(currentLoadCase)%P%maskLogical .and. transpose(loadCases(currentLoadCase)%P%maskLogical) .and. & if (any(loadCases(currentLoadCase)%P%maskLogical .and. transpose(loadCases(currentLoadCase)%P%maskLogical) .and. &
reshape([ .false.,.true.,.true.,.true.,.false.,.true.,.true.,.true.,.false.],[ 3,3]))) & reshape([ .false.,.true.,.true.,.true.,.false.,.true.,.true.,.true.,.false.],[ 3,3]))) &
@ -290,25 +290,29 @@ program DAMASK_spectral_Driver
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! write header of output file ! write header of output file
if (appendToOutFile) then if (appendToOutFile) then
open(538,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.spectralOut',& open(newunit=resUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
form='UNFORMATTED', position='APPEND', status='OLD') '.spectralOut',form='UNFORMATTED', position='APPEND', status='OLD')
open(newunit=statUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
'.sta',form='FORMATTED', position='APPEND', status='OLD')
else else
open(538,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//'.spectralOut',& open(newunit=resUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
form='UNFORMATTED',status='REPLACE') '.spectralOut',form='UNFORMATTED',status='REPLACE')
write(538) 'load', trim(loadCaseFile) open(newunit=statUnit,file=trim(getSolverWorkingDirectoryName())//trim(getSolverJobName())//&
write(538) 'workingdir', trim(getSolverWorkingDirectoryName()) '.sta',form='FORMATTED',status='REPLACE')
write(538) 'geometry', trim(geometryFile) write(resUnit) 'load', trim(loadCaseFile)
write(538) 'resolution', res write(resUnit) 'workingdir', trim(getSolverWorkingDirectoryName())
write(538) 'dimension', geomdim write(resUnit) 'geometry', trim(geometryFile)
write(538) 'materialpoint_sizeResults', materialpoint_sizeResults write(resUnit) 'resolution', res
write(538) 'loadcases', size(loadCases) write(resUnit) 'dimension', geomdim
write(538) 'frequencies', loadCases%outputfrequency ! one entry per currentLoadCase write(resUnit) 'materialpoint_sizeResults', materialpoint_sizeResults
write(538) 'times', loadCases%time ! one entry per currentLoadCase write(resUnit) 'loadcases', size(loadCases)
write(538) 'logscales', loadCases%logscale write(resUnit) 'frequencies', loadCases%outputfrequency ! one entry per currentLoadCase
write(538) 'increments', loadCases%incs ! one entry per currentLoadCase write(resUnit) 'times', loadCases%time ! one entry per currentLoadCase
write(538) 'startingIncrement', restartInc - 1_pInt ! start with writing out the previous inc write(resUnit) 'logscales', loadCases%logscale
write(538) 'eoh' ! end of header write(resUnit) 'increments', loadCases%incs ! one entry per currentLoadCase
write(538) materialpoint_results(1_pInt:materialpoint_sizeResults,1,1_pInt:mesh_NcpElems) ! initial (non-deformed or read-in) results write(resUnit) 'startingIncrement', restartInc - 1_pInt ! start with writing out the previous inc
write(resUnit) 'eoh' ! end of header
write(resUnit) materialpoint_results(1_pInt:materialpoint_sizeResults,1,1_pInt:mesh_NcpElems) ! initial (non-deformed or read-in) results
if (debugGeneral) write(6,'(a)') 'Header of result file written out' if (debugGeneral) write(6,'(a)') 'Header of result file written out'
endif endif
@ -408,11 +412,13 @@ program DAMASK_spectral_Driver
elseif (solres%termIll) then ! material point model cannot find a solution elseif (solres%termIll) then ! material point model cannot find a solution
call IO_error(850_pInt) call IO_error(850_pInt)
else else
guessmode=1.0_pReal ! start guessing after first accepted (not converged) (sub)inc guessmode = 1.0_pReal ! start guessing after first accepted (not converged) (sub)inc
endif endif
else else
guessmode = 1.0_pReal ! start guessing after first converged (sub)inc guessmode = 1.0_pReal ! start guessing after first converged (sub)inc
endif endif
if(guessmode == 1.0_pReal) &
write(statUnit,*) inc, time, cutBackLevel, solres%converged, solres%iterationsNeeded
enddo subIncLooping enddo subIncLooping
cutBackLevel = max(0_pInt, cutBackLevel - 1_pInt) ! try half subincs next inc cutBackLevel = max(0_pInt, cutBackLevel - 1_pInt) ! try half subincs next inc
if(solres%converged) then if(solres%converged) then
@ -427,7 +433,7 @@ program DAMASK_spectral_Driver
if (mod(inc,loadCases(currentLoadCase)%outputFrequency) == 0_pInt) then ! at output frequency if (mod(inc,loadCases(currentLoadCase)%outputFrequency) == 0_pInt) then ! at output frequency
write(6,'(1/,a)') '... writing results to file ......................................' write(6,'(1/,a)') '... writing results to file ......................................'
write(538) materialpoint_results ! write result to file write(resUnit) materialpoint_results ! write result to file
endif endif
else !just time forwarding else !just time forwarding
time = time + timeinc time = time + timeinc
@ -456,7 +462,7 @@ program DAMASK_spectral_Driver
real(convergedCounter, pReal)/& real(convergedCounter, pReal)/&
real(notConvergedCounter + convergedCounter,pReal)*100.0_pReal, & real(notConvergedCounter + convergedCounter,pReal)*100.0_pReal, &
' %) increments converged!' ' %) increments converged!'
close(538) close(resUnit)
if (notConvergedCounter > 0_pInt) call quit(3_pInt) if (notConvergedCounter > 0_pInt) call quit(3_pInt)
call quit(0_pInt) call quit(0_pInt)

10
code/DAMASK_spectral_SolverAL.f90 → code/DAMASK_spectral_solverAL.f90
View File

@ -71,7 +71,7 @@ module DAMASK_spectral_SolverAL
real(pReal), private :: err_stress, err_f, err_p real(pReal), private :: err_stress, err_f, err_p
logical, private :: ForwardData logical, private :: ForwardData
real(pReal), private, dimension(3,3) :: mask_stress = 0.0_pReal real(pReal), private, dimension(3,3) :: mask_stress = 0.0_pReal
integer(pInt) :: reportIter = 0_pInt
contains contains
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -343,8 +343,10 @@ else
terminallyIll = .false. terminallyIll = .false.
if (reason > 0 ) AL_solution%converged = .true. if (reason > 0 ) then
AL_solution%converged = .true.
AL_solution%iterationsNeeded = reportIter
endif
end function AL_solution end function AL_solution
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -373,7 +375,7 @@ else
implicit none implicit none
integer(pInt) :: i,j,k integer(pInt) :: i,j,k
integer(pInt), save :: callNo = 3_pInt, reportIter = 0_pInt integer(pInt), save :: callNo = 3_pInt
real(pReal), dimension(3,3) :: temp33_Real real(pReal), dimension(3,3) :: temp33_Real
logical :: report logical :: report

14
code/DAMASK_spectral_SolverBasic.f90 → code/DAMASK_spectral_solverBasic.f90
View File

@ -235,14 +235,15 @@ else
F_lastInc = F F_lastInc = F
endif endif
F_aim = F_aim + f_aimDot * timeinc F_aim = F_aim + f_aimDot * timeinc
F = Utilities_forwardField(timeinc,F_aim,F_lastInc,Fdot) !I thin F aim should be rotated here F = Utilities_forwardField(timeinc,F_aim,F_lastInc,Fdot) !I think F aim should be rotated here
print*, 'F', sum(sum(sum(F,dim=5),dim=4),dim=3)*wgt
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! update stiffness (and gamma operator) ! update stiffness (and gamma operator)
S = Utilities_maskedCompliance(rotation_BC,P_BC%maskLogical,C) S = Utilities_maskedCompliance(rotation_BC,P_BC%maskLogical,C)
if (update_gamma) call Utilities_updateGamma(C) if (update_gamma) call Utilities_updateGamma(C)
ForwardData = .True.
iter = 0_pInt iter = 0_pInt
convergenceLoop: do while(iter < itmax) convergenceLoop: do while(iter < itmax)
@ -250,7 +251,7 @@ print*, 'F', sum(sum(sum(F,dim=5),dim=4),dim=3)*wgt
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! report begin of new iteration ! report begin of new iteration
write(6,'(/,a,3(a,'//IO_intOut(itmax)//'))') trim(incInfo), & write(6,'(/,a,3(a,'//IO_intOut(itmax)//'))') trim(incInfo), &
' @ Iter. ', itmin, '<',iter, '≤', itmax ' @ Iter. ', itmin, '≤',iter, '≤', itmax
write(6,'(a,/,3(3(f12.7,1x)/))',advance='no') 'deformation gradient aim =', & write(6,'(a,/,3(3(f12.7,1x)/))',advance='no') 'deformation gradient aim =', &
math_transpose33(F_aim) math_transpose33(F_aim)
F_aim_lab_lastIter = math_rotate_backward33(F_aim,rotation_BC) F_aim_lab_lastIter = math_rotate_backward33(F_aim,rotation_BC)
@ -282,7 +283,10 @@ print*, 'F', sum(sum(sum(F,dim=5),dim=4),dim=3)*wgt
F = F - reshape(field_real(1:res(1),1:res(2),1:res(3),1:3,1:3),shape(F),order=[3,4,5,1,2]) ! F(x)^(n+1) = F(x)^(n) + correction; *wgt: correcting for missing normalization F = F - reshape(field_real(1:res(1),1:res(2),1:res(3),1:3,1:3),shape(F),order=[3,4,5,1,2]) ! F(x)^(n+1) = F(x)^(n) + correction; *wgt: correcting for missing normalization
basic_solution%converged = basic_Converged(err_div,P_av,err_stress,P_av) basic_solution%converged = basic_Converged(err_div,P_av,err_stress,P_av)
write(6,'(/,a)') '==========================================================================' write(6,'(/,a)') '=========================================================================='
if ((basic_solution%converged .and. iter > itmin) .or. basic_solution%termIll) exit if ((basic_solution%converged .and. iter >= itmin) .or. basic_solution%termIll) then
basic_solution%iterationsNeeded = iter
exit
endif
enddo convergenceLoop enddo convergenceLoop
end function basic_solution end function basic_solution

7
code/DAMASK_spectral_SolverBasicPETSC.f90 → code/DAMASK_spectral_solverBasicPETSc.f90
View File

@ -116,7 +116,7 @@ subroutine BasicPETSC_init()
call Utilities_init() call Utilities_init()
write(6,'(/,a)') ' <<<+- DAMASK_spectral_solverBasicPETSC init -+>>>' write(6,'(/,a)') ' <<<+- DAMASK_spectral_solverBasicPETSc init -+>>>'
write(6,'(a)') ' $Id: DAMASK_spectral_SolverBasicPETSC.f90 1654 2012-08-03 09:25:48Z MPIE\m.diehl $' write(6,'(a)') ' $Id: DAMASK_spectral_SolverBasicPETSC.f90 1654 2012-08-03 09:25:48Z MPIE\m.diehl $'
#include "compilation_info.f90" #include "compilation_info.f90"
write(6,'(a)') '' write(6,'(a)') ''
@ -307,7 +307,9 @@ else
basicPETSc_solution%termIll = terminallyIll basicPETSc_solution%termIll = terminallyIll
terminallyIll = .false. terminallyIll = .false.
BasicPETSC_solution%converged =.false. BasicPETSC_solution%converged =.false.
if (reason > 0 ) BasicPETSC_solution%converged = .true. if (reason > 0 ) then
BasicPETSC_solution%converged = .true.
endif
end function BasicPETSC_solution end function BasicPETSC_solution
@ -433,7 +435,6 @@ else
' (',err_div/pAvgDivL2,' N/m³)' ' (',err_div/pAvgDivL2,' N/m³)'
write(6,'(a,f6.2,a,es11.4,a)') 'error stress = ', err_stress/min(maxval(abs(P_av))*err_stress_tolrel,err_stress_tolabs), & write(6,'(a,f6.2,a,es11.4,a)') 'error stress = ', err_stress/min(maxval(abs(P_av))*err_stress_tolrel,err_stress_tolabs), &
' (',err_stress,' Pa)' ' (',err_stress,' Pa)'
return
end subroutine BasicPETSC_converged end subroutine BasicPETSC_converged

46
code/DAMASK_spectral_Utilities.f90 → code/DAMASK_spectral_utilities.f90
View File

@ -23,12 +23,11 @@ module DAMASK_spectral_Utilities
IO_error IO_error
implicit none implicit none
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! variables storing information for spectral method and FFTW ! variables storing information for spectral method and FFTW
type(C_PTR), private :: plan_forward, plan_backward ! plans for fftw type(C_PTR), private :: plan_forward, plan_backward ! plans for fftw
real(pReal), private, dimension(:,:,:,:,:,:,:), allocatable :: gamma_hat ! gamma operator (field) for spectral method real(pReal), private, dimension(:,:,:,:,:,:,:), allocatable :: gamma_hat ! gamma operator (field) for spectral method
real(pReal), private, dimension(:,:,:,:), allocatable :: xi ! wave vector field for divergence and for gamma operator real(pReal), private, dimension(:,:,:,:), allocatable :: xi ! wave vector field for divergence and for gamma operator
complex(pReal),private, dimension(:,:,:,:,:), pointer :: field_fourier complex(pReal),private, dimension(:,:,:,:,:), pointer :: field_fourier
real(pReal), private, dimension(3,3,3,3) :: C_ref real(pReal), private, dimension(3,3,3,3) :: C_ref
@ -54,9 +53,10 @@ module DAMASK_spectral_Utilities
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! derived types ! derived types
type tSolutionState type tSolutionState
logical :: converged = .true. logical :: converged = .true.
logical :: regrid = .false. logical :: regrid = .false.
logical :: termIll = .false. logical :: termIll = .false.
integer(pInt) :: iterationsNeeded = 0_pInt
end type tSolutionState end type tSolutionState
type tBoundaryCondition type tBoundaryCondition
@ -105,7 +105,7 @@ subroutine Utilities_init()
write(6,'(a)') '' write(6,'(a)') ''
write(6,'(a)') ' <<<+- DAMASK_Utilities init -+>>>' write(6,'(a)') ' <<<+- DAMASK_spectral_utilities init -+>>>'
write(6,'(a)') ' $Id$' write(6,'(a)') ' $Id$'
#include "compilation_info.f90" #include "compilation_info.f90"
write(6,'(a)') '' write(6,'(a)') ''
@ -454,14 +454,15 @@ end function Utilities_divergenceRMS
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates mask compliance !> @brief calculates mask compliance
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
function Utilities_maskedCompliance(rot_BC,mask_stressVector,C) function Utilities_maskedCompliance(rot_BC,mask_stress,C)
implicit none implicit none
real(pReal), dimension(3,3,3,3) :: Utilities_maskedCompliance real(pReal), dimension(3,3,3,3) :: Utilities_maskedCompliance
real(pReal), dimension(3,3,3,3), intent(in) :: C real(pReal), dimension(3,3,3,3), intent(in) :: C
integer(pInt) :: j, k, m, n integer(pInt) :: j, k, m, n
real(pReal), dimension(3,3), intent(in) :: rot_BC real(pReal), dimension(3,3), intent(in) :: rot_BC
logical, dimension(9), intent(in) :: mask_stressVector logical, dimension(3,3), intent(in) :: mask_stress
logical, dimension(9) :: mask_stressVector
real(pReal), dimension(3,3,3,3) :: C_lastInc real(pReal), dimension(3,3,3,3) :: C_lastInc
real(pReal), dimension(9,9) :: temp99_Real real(pReal), dimension(9,9) :: temp99_Real
integer(pInt) :: size_reduced = 0_pInt integer(pInt) :: size_reduced = 0_pInt
@ -469,13 +470,16 @@ function Utilities_maskedCompliance(rot_BC,mask_stressVector,C)
logical :: errmatinv logical :: errmatinv
character(len=1024):: formatString character(len=1024):: formatString
size_reduced = count(mask_stressVector) mask_stressVector = reshape(transpose(mask_stress), [9])
size_reduced = int(count(mask_stressVector), pInt)
if(size_reduced > 0_pInt )then if(size_reduced > 0_pInt )then
allocate (c_reduced(size_reduced,size_reduced), source =0.0_pReal) allocate (c_reduced(size_reduced,size_reduced), source =0.0_pReal)
allocate (s_reduced(size_reduced,size_reduced), source =0.0_pReal) allocate (s_reduced(size_reduced,size_reduced), source =0.0_pReal)
allocate (sTimesC(size_reduced,size_reduced), source =0.0_pReal) allocate (sTimesC(size_reduced,size_reduced), source =0.0_pReal)
C_lastInc = math_rotate_forward3333(C*0.75_pReal+0.05_pReal*C_ref,rot_BC) ! calculate stiffness from former inc C_lastInc = math_rotate_forward3333(C,rot_BC) ! calculate stiffness from former inc
print*,'C'
print*,C_lastInc
temp99_Real = math_Plain3333to99(C_lastInc) temp99_Real = math_Plain3333to99(C_lastInc)
k = 0_pInt ! build reduced stiffness k = 0_pInt ! build reduced stiffness
do n = 1_pInt,9_pInt do n = 1_pInt,9_pInt
@ -521,7 +525,8 @@ function Utilities_maskedCompliance(rot_BC,mask_stressVector,C)
temp99_real = 0.0_pReal temp99_real = 0.0_pReal
endif endif
Utilities_maskedCompliance = math_Plain99to3333(temp99_Real) Utilities_maskedCompliance = math_Plain99to3333(temp99_Real)
print*,'masked S'
print*,Utilities_maskedCompliance
end function Utilities_maskedCompliance end function Utilities_maskedCompliance
subroutine Utilities_constitutiveResponse(coordinates,F_lastInc,F,temperature,timeinc,& subroutine Utilities_constitutiveResponse(coordinates,F_lastInc,F,temperature,timeinc,&
@ -539,15 +544,16 @@ subroutine Utilities_constitutiveResponse(coordinates,F_lastInc,F,temperature,ti
real(pReal), dimension(3,3,res(1),res(2),res(3)), intent(in) :: F,F_lastInc real(pReal), dimension(3,3,res(1),res(2),res(3)), intent(in) :: F,F_lastInc
real(pReal), dimension(3,3,res(1),res(2),res(3)) :: P real(pReal), dimension(3,3,res(1),res(2),res(3)) :: P
real(pReal) :: timeinc real(pReal),intent(in) :: timeinc
logical, intent(in) :: forwardData logical, intent(in) :: forwardData
integer(pInt) :: i, j, k, ielem integer(pInt) :: i, j, k, ielem
integer(pInt) :: calcMode, collectMode integer(pInt) :: calcMode, collectMode
real(pReal), dimension(3,3,3,3) :: dPdF, C real(pReal), dimension(3,3,3,3) :: dPdF
real(pReal), dimension(3,3,3,3),intent(out) :: C
real(pReal), dimension(6) :: sigma ! cauchy stress real(pReal), dimension(6) :: sigma ! cauchy stress
real(pReal), dimension(6,6) :: dsde real(pReal), dimension(6,6) :: dsde
real(pReal), dimension(3,3), intent(in) :: rotation_BC real(pReal), dimension(3,3), intent(in) :: rotation_BC
real(pReal), dimension(3,3) :: P_av real(pReal), dimension(3,3),intent(out) :: P_av
write(6,'(/,a,/)') '... evaluating constitutive response ......................................' write(6,'(/,a,/)') '... evaluating constitutive response ......................................'
if (forwardData) then if (forwardData) then
@ -561,8 +567,12 @@ subroutine Utilities_constitutiveResponse(coordinates,F_lastInc,F,temperature,ti
calcMode = 2_pInt calcMode = 2_pInt
collectMode = 5_pInt collectMode = 5_pInt
endif endif
print*, 'collect mode', collectMode
print*, 'calc mode', calcMode if (DebugGeneral) then
write(6,*) 'collect mode', collectMode
write(6,*) 'calc mode', calcMode
endif
ielem = 0_pInt ielem = 0_pInt
do k = 1_pInt, res(3); do j = 1_pInt, res(2); do i = 1_pInt, res(1) do k = 1_pInt, res(3); do j = 1_pInt, res(2); do i = 1_pInt, res(1)
ielem = ielem + 1_pInt ielem = ielem + 1_pInt

36
code/Makefile
View File

@ -127,6 +127,7 @@ endif
ifneq "$(FASTBUILD)" "YES" ifneq "$(FASTBUILD)" "YES"
STANDARD_CHECK_ifort ?=-stand f08 -standard-semantics -warn stderrors STANDARD_CHECK_ifort ?=-stand f08 -standard-semantics -warn stderrors
STANDARD_CHECK_gfortran ?=-std=f2008 -fall-intrinsics STANDARD_CHECK_gfortran ?=-std=f2008 -fall-intrinsics
#-std=2008ts for newer gfortran
endif endif
#-fall-intrinsics: all intrinsic procedures (including the GNU-specific extensions) are accepted. This can be useful with -std=f95 to force standard-compliance #-fall-intrinsics: all intrinsic procedures (including the GNU-specific extensions) are accepted. This can be useful with -std=f95 to force standard-compliance
# but get access to the full range of intrinsics available with gfortran. As a consequence, -Wintrinsics-std will be ignored and no user-defined # but get access to the full range of intrinsics available with gfortran. As a consequence, -Wintrinsics-std will be ignored and no user-defined
@ -142,6 +143,7 @@ OPTIMIZATION_AGGRESSIVE_gfortran :=-O3 $(PORTABLE_SWITCH) -ffast-math -funroll-l
COMPILE_OPTIONS_ifort :=-fpp\ COMPILE_OPTIONS_ifort :=-fpp\
-ftz\
-assume byterecl -assume byterecl
ifneq "$(FASTBUILD)" "YES" ifneq "$(FASTBUILD)" "YES"
COMPILE_OPTIONS_ifort :=$(COMPILE_OPTIONS_ifort)\ COMPILE_OPTIONS_ifort :=$(COMPILE_OPTIONS_ifort)\
@ -159,6 +161,7 @@ COMPILE_OPTIONS_ifort :=$(COMPILE_OPTIONS_ifort)\
endif endif
#-fpp: preprocessor #-fpp: preprocessor
#-ftz: flush unterflow to zero, automatically set if O<0,1,2,3> >0
#-fimplicit-none: assume "implicit-none" even if not present in source #-fimplicit-none: assume "implicit-none" even if not present in source
#-diag-disable: disables warnings, where #-diag-disable: disables warnings, where
# warning ID 5268: the text exceeds right hand column allowed on the line (we have only comments there) # warning ID 5268: the text exceeds right hand column allowed on the line (we have only comments there)
@ -194,6 +197,7 @@ endif
# pointers: Checking for certain disassociated or uninitialized pointers or unallocated allocatable objects. # pointers: Checking for certain disassociated or uninitialized pointers or unallocated allocatable objects.
# uninit: Checking for uninitialized variables. # uninit: Checking for uninitialized variables.
#-heap-arrays: should not be done for OpenMP, but set "ulimit -s unlimited" on shell. Probably it helps also to unlimit other limits #-heap-arrays: should not be done for OpenMP, but set "ulimit -s unlimited" on shell. Probably it helps also to unlimit other limits
#-fpe-all0 capture all floating-point exceptions, sets -ftz automatically
################################################################################################### ###################################################################################################
COMPILE_OPTIONS_gfortran :=-xf95-cpp-input COMPILE_OPTIONS_gfortran :=-xf95-cpp-input
@ -258,6 +262,12 @@ endif
# #
#OPTIONS FOR DEGUBBING DURING RUNTIME #OPTIONS FOR DEGUBBING DURING RUNTIME
#-fcheck-bounds: check if an array index is too small (<1) or too large! #-fcheck-bounds: check if an array index is too small (<1) or too large!
#-ffpe-trap=invalid,\ stop execution if floating point exception is detected (NaN is silent)
# zero,\
# overflow,\
# underflow,\
# precision,\
# denormal
################################################################################################## ##################################################################################################
PRECISION_ifort :=-real-size 64 -integer-size 32 -DFLOAT=8 -DINT=4 PRECISION_ifort :=-real-size 64 -integer-size 32 -DFLOAT=8 -DINT=4
@ -281,29 +291,29 @@ COMPILED_FILES = prec.o DAMASK_spectral_interface.o IO.o numerics.o debug.o math
ifeq "$(SOLVER)" "NEW" ifeq "$(SOLVER)" "NEW"
ifdef PETSC_DIR ifdef PETSC_DIR
PETSC_FILES = DAMASK_spectral_SolverAL.o DAMASK_spectral_SolverBasicPETSC.o PETSC_FILES = DAMASK_spectral_solverAL.o DAMASK_spectral_solverBasicPETSc.o
endif endif
COMPILED_FILES += DAMASK_spectral_Utilities.o DAMASK_spectral_SolverBasic.o $(PETSC_FILES) COMPILED_FILES += DAMASK_spectral_utilities.o DAMASK_spectral_solverBasic.o $(PETSC_FILES)
DAMASK_spectral.exe: DAMASK_spectral_Driver.o DAMASK_spectral.exe: DAMASK_spectral_driver.o
$(PREFIX) $(COMPILERNAME) $(OPENMP_FLAG_$(F90)) $(OPTIMIZATION_$(MAXOPTI)_$(F90)) $(STANDARD_CHECK_$(F90)) \ $(PREFIX) $(COMPILERNAME) $(OPENMP_FLAG_$(F90)) $(OPTIMIZATION_$(MAXOPTI)_$(F90)) $(STANDARD_CHECK_$(F90)) \
-o DAMASK_spectral.exe DAMASK_spectral_Driver.o \ -o DAMASK_spectral.exe DAMASK_spectral_driver.o \
$(COMPILED_FILES) $(LIB_DIRS) $(LIBRARIES) $(SUFFIX) $(COMPILED_FILES) $(LIB_DIRS) $(LIBRARIES) $(SUFFIX)
DAMASK_spectral_Driver.o: DAMASK_spectral_Driver.f90 DAMASK_spectral_SolverBasic.o $(PETSC_FILES) DAMASK_spectral_driver.o: DAMASK_spectral_driver.f90 DAMASK_spectral_solverBasic.o $(PETSC_FILES)
$(PREFIX) $(COMPILERNAME) $(COMPILE_MAXOPTI) -c DAMASK_spectral_Driver.f90 $(SUFFIX) $(PREFIX) $(COMPILERNAME) $(COMPILE_MAXOPTI) -c DAMASK_spectral_driver.f90 $(SUFFIX)
DAMASK_spectral_SolverAL.o: DAMASK_spectral_SolverAL.f90\ DAMASK_spectral_solverAL.o: DAMASK_spectral_solverAL.f90\
DAMASK_spectral_Utilities.o DAMASK_spectral_utilities.o
DAMASK_spectral_SolverBasic.o: DAMASK_spectral_SolverBasic.f90\ DAMASK_spectral_solverBasic.o: DAMASK_spectral_solverBasic.f90\
DAMASK_spectral_Utilities.o DAMASK_spectral_utilities.o
DAMASK_spectral_SolverBasicPETSC.o: DAMASK_spectral_SolverBasicPETSC.f90\ DAMASK_spectral_solverBasicPETSc.o: DAMASK_spectral_solverBasicPETSc.f90\
DAMASK_spectral_Utilities.o DAMASK_spectral_utilities.o
DAMASK_spectral_Utilities.o: DAMASK_spectral_Utilities.f90\ DAMASK_spectral_utilities.o: DAMASK_spectral_utilities.f90\
CPFEM.o CPFEM.o
else else
DAMASK_spectral.exe: DAMASK_spectral.o DAMASK_spectral.exe: DAMASK_spectral.o