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

removed deadlock when calling IO_warning(33) 

polished IO_error and IO_warning
checked OpenMP critical statements in DAMASK_spectral.f90
corrected writing of headed in DAMASK_spectral.f90
This commit is contained in:
Martin Diehl 2011-11-02 14:38:42 +00:00
parent 447e69019c
commit 680ba9082f
5 changed files with 108 additions and 78 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

79
code/DAMASK_spectral.f90
View File

@ -85,6 +85,7 @@ program DAMASK_spectral
bc_velGradApplied ! decide wether velocity gradient or fdot is given bc_velGradApplied ! decide wether velocity gradient or fdot is given
logical, dimension(:,:,:,:), allocatable :: bc_mask ! mask of boundary conditions logical, dimension(:,:,:,:), allocatable :: bc_mask ! mask of boundary conditions
logical, dimension(:,:,:), allocatable :: bc_maskvector ! linear mask of boundary conditions logical, dimension(:,:,:), allocatable :: bc_maskvector ! linear mask of boundary conditions
character(len=3) :: loadcase_string
! variables storing information from geom file ! variables storing information from geom file
real(pReal) :: wgt real(pReal) :: wgt
@ -152,7 +153,7 @@ program DAMASK_spectral
resolution = 1_pInt resolution = 1_pInt
geomdimension = 0.0_pReal geomdimension = 0.0_pReal
if (command_argument_count() /= 4) call IO_error(102) ! check for correct number of given arguments if (command_argument_count() /= 4) call IO_error(error_ID=102) ! check for correct number of given arguments
! Reading the loadcase file and allocate variables ! Reading the loadcase file and allocate variables
path = getLoadcaseName() path = getLoadcaseName()
@ -162,7 +163,7 @@ program DAMASK_spectral
print '(a,a)', 'Solver Job Name: ',trim(getSolverJobName()) print '(a,a)', 'Solver Job Name: ',trim(getSolverJobName())
print '(a)', '******************************************************' print '(a)', '******************************************************'
!$OMP END CRITICAL (write2out) !$OMP END CRITICAL (write2out)
if (.not. IO_open_file(myUnit,path)) call IO_error(30,ext_msg = trim(path)) if (.not. IO_open_file(myUnit,path)) call IO_error(error_ID=30,ext_msg = trim(path))
rewind(myUnit) rewind(myUnit)
do do
@ -185,7 +186,7 @@ program DAMASK_spectral
100 N_Loadcases = N_n 100 N_Loadcases = N_n
if ((N_l + N_Fdot /= N_n) .or. (N_n /= N_t)) & ! sanity check if ((N_l + N_Fdot /= N_n) .or. (N_n /= N_t)) & ! sanity check
call IO_error(37,ext_msg = trim(path)) ! error message for incomplete loadcase call IO_error(error_ID=37,ext_msg = trim(path)) ! error message for incomplete loadcase
allocate (bc_deformation(3,3,N_Loadcases)); bc_deformation = 0.0_pReal allocate (bc_deformation(3,3,N_Loadcases)); bc_deformation = 0.0_pReal
allocate (bc_stress(3,3,N_Loadcases)); bc_stress = 0.0_pReal allocate (bc_stress(3,3,N_Loadcases)); bc_stress = 0.0_pReal
@ -270,7 +271,7 @@ program DAMASK_spectral
path = getModelName() path = getModelName()
if (.not. IO_open_file(myUnit,trim(path)//InputFileExtension))& if (.not. IO_open_file(myUnit,trim(path)//InputFileExtension))&
call IO_error(101,ext_msg = trim(path)//InputFileExtension) call IO_error(error_ID=101,ext_msg = trim(path)//InputFileExtension)
rewind(myUnit) rewind(myUnit)
read(myUnit,'(a1024)') line read(myUnit,'(a1024)') line
@ -279,7 +280,7 @@ program DAMASK_spectral
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,posGeom,1) + 1_pInt headerLength = IO_intValue(line,posGeom,1) + 1_pInt
else else
call IO_error(42) call IO_error(error_ID=42)
endif endif
rewind(myUnit) rewind(myUnit)
@ -319,11 +320,11 @@ program DAMASK_spectral
end select end select
enddo enddo
close(myUnit) close(myUnit)
if (.not.(gotDimension .and. gotHomogenization .and. gotResolution)) call IO_error(45) if (.not.(gotDimension .and. gotHomogenization .and. gotResolution)) call IO_error(error_ID=45)
if(mod(resolution(1),2_pInt)/=0_pInt .or.& if(mod(resolution(1),2_pInt)/=0_pInt .or.&
mod(resolution(2),2_pInt)/=0_pInt .or.& mod(resolution(2),2_pInt)/=0_pInt .or.&
(mod(resolution(3),2_pInt)/=0_pInt .and. resolution(3)/= 1_pInt)) call IO_error(103) (mod(resolution(3),2_pInt)/=0_pInt .and. resolution(3)/= 1_pInt)) call IO_error(error_ID=103)
allocate (defgrad ( resolution(1),resolution(2),resolution(3),3,3)); defgrad = 0.0_pReal allocate (defgrad ( resolution(1),resolution(2),resolution(3),3,3)); defgrad = 0.0_pReal
allocate (defgradold ( resolution(1),resolution(2),resolution(3),3,3)); defgradold = 0.0_pReal allocate (defgradold ( resolution(1),resolution(2),resolution(3),3,3)); defgradold = 0.0_pReal
@ -350,52 +351,68 @@ program DAMASK_spectral
print '(a,L)','spectralPictureMode: ',spectralPictureMode print '(a,L)','spectralPictureMode: ',spectralPictureMode
print '(a)', '******************************************************' print '(a)', '******************************************************'
print '(a,a)','Loadcase File Name: ',trim(getLoadcaseName()) print '(a,a)','Loadcase File Name: ',trim(getLoadcaseName())
!$OMP END CRITICAL (write2out)
if (bc_followFormerTrajectory(1)) then if (bc_followFormerTrajectory(1)) then
call IO_warning(33) ! cannot guess along trajectory for first step of first loadcase call IO_warning(warning_ID=33_pInt) ! cannot guess along trajectory for first step of first loadcase
bc_followFormerTrajectory(1) = .false. bc_followFormerTrajectory(1) = .false.
endif endif
! consistency checks and output of loadcase ! consistency checks and output of loadcase
do loadcase = 1, N_Loadcases do loadcase = 1, N_Loadcases
!$OMP CRITICAL (write2out)
print '(a)', '------------------------------------------------------' print '(a)', '------------------------------------------------------'
print '(a,i5)', 'Loadcase: ', loadcase print '(a,i5)', 'Loadcase: ', loadcase
write (loadcase_string, '(i3)' ) loadcase
if (.not. bc_followFormerTrajectory(loadcase)) & if (.not. bc_followFormerTrajectory(loadcase)) &
print '(a)', 'Drop Guessing Along Trajectory' print '(a)', 'Drop Guessing Along Trajectory'
!$OMP END CRITICAL (write2out)
if (any(bc_mask(:,:,1,loadcase) .eqv. bc_mask(1:3,1:3,2,loadcase)))& ! exclusive or masking only if (any(bc_mask(:,:,1,loadcase) .eqv. bc_mask(1:3,1:3,2,loadcase)))& ! exclusive or masking only
call IO_error(31,loadcase) call IO_error(error_ID=31,ext_msg=loadcase_string)
if (any(bc_mask(1:3,1:3,2,loadcase).and.transpose(bc_mask(1:3,1:3,2,loadcase)).and.& !checking if no rotation is allowed by stress BC if (any(bc_mask(1:3,1:3,2,loadcase).and.transpose(bc_mask(1:3,1:3,2,loadcase)).and.& !checking if no rotation is allowed by stress BC
reshape((/.false.,.true.,.true.,.true.,.false.,.true.,.true.,.true.,.false./),(/3,3/))))& reshape((/.false.,.true.,.true.,.true.,.false.,.true.,.true.,.true.,.false./),(/3,3/))))&
call IO_error(38,loadcase) call IO_error(error_ID=38,ext_msg=loadcase_string)
if (bc_velGradApplied(loadcase)) then if (bc_velGradApplied(loadcase)) then
do j = 1, 3 do j = 1, 3
if (any(bc_mask(j,1:3,1,loadcase) .eqv. .true.) .and.& if (any(bc_mask(j,1:3,1,loadcase) .eqv. .true.) .and.&
any(bc_mask(j,1:3,1,loadcase) .eqv. .false.)) call IO_error(32,loadcase) ! each line should be either fully or not at all defined any(bc_mask(j,1:3,1,loadcase) .eqv. .false.)) call IO_error(error_ID=32,ext_msg=loadcase_string) ! each line should be either fully or not at all defined
enddo enddo
!$OMP CRITICAL (write2out)
print '(a,/,3(3(f12.6,x)/))','Velocity Gradient:', merge(math_transpose3x3(bc_deformation(1:3,1:3,loadcase)),& print '(a,/,3(3(f12.6,x)/))','Velocity Gradient:', merge(math_transpose3x3(bc_deformation(1:3,1:3,loadcase)),&
reshape(spread(DAMASK_NaN,1,9),(/3,3/)),& reshape(spread(DAMASK_NaN,1,9),(/3,3/)),&
transpose(bc_mask(1:3,1:3,1,loadcase))) transpose(bc_mask(1:3,1:3,1,loadcase)))
!$OMP END CRITICAL (write2out)
else else
!$OMP CRITICAL (write2out)
print '(a,/,3(3(f12.6,x)/))','Change of Deformation Gradient:', merge(math_transpose3x3(bc_deformation(1:3,1:3,loadcase)),& print '(a,/,3(3(f12.6,x)/))','Change of Deformation Gradient:', merge(math_transpose3x3(bc_deformation(1:3,1:3,loadcase)),&
reshape(spread(DAMASK_NaN,1,9),(/3,3/)),& reshape(spread(DAMASK_NaN,1,9),(/3,3/)),&
transpose(bc_mask(1:3,1:3,1,loadcase))) transpose(bc_mask(1:3,1:3,1,loadcase)))
!$OMP END CRITICAL (write2out)
endif endif
!$OMP CRITICAL (write2out)
print '(a,/,3(3(f12.6,x)/))','Stress Boundary Condition/MPa:',merge(math_transpose3x3(bc_stress(1:3,1:3,loadcase)),& print '(a,/,3(3(f12.6,x)/))','Stress Boundary Condition/MPa:',merge(math_transpose3x3(bc_stress(1:3,1:3,loadcase)),&
reshape(spread(DAMASK_NaN,1,9),(/3,3/)),& reshape(spread(DAMASK_NaN,1,9),(/3,3/)),&
transpose(bc_mask(:,:,2,loadcase)))*1e-6 transpose(bc_mask(:,:,2,loadcase)))*1e-6
!$OMP END CRITICAL (write2out)
if (any(abs(math_mul33x33(bc_rotation(1:3,1:3,loadcase),math_transpose3x3(bc_rotation(1:3,1:3,loadcase)))-math_I3)& if (any(abs(math_mul33x33(bc_rotation(1:3,1:3,loadcase),math_transpose3x3(bc_rotation(1:3,1:3,loadcase)))-math_I3)&
>reshape(spread(rotation_tol,1,9),(/3,3/)))& >reshape(spread(rotation_tol,1,9),(/3,3/)))&
.or. abs(math_det3x3(bc_rotation(1:3,1:3,loadcase)))>1.0_pReal + rotation_tol) call IO_error(46,loadcase) .or. abs(math_det3x3(bc_rotation(1:3,1:3,loadcase)))>1.0_pReal + rotation_tol) call IO_error(error_ID=46,ext_msg=loadcase_string)
!$OMP CRITICAL (write2out)
if (any(bc_rotation(1:3,1:3,loadcase)/=math_I3)) & if (any(bc_rotation(1:3,1:3,loadcase)/=math_I3)) &
print '(a,/,3(3(f12.6,x)/))','Rotation of BCs:',math_transpose3x3(bc_rotation(1:3,1:3,loadcase)) print '(a,/,3(3(f12.6,x)/))','Rotation of BCs:',math_transpose3x3(bc_rotation(1:3,1:3,loadcase))
if (bc_timeIncrement(loadcase) < 0.0_pReal) call IO_error(34,loadcase) ! negative time increment !$OMP END CRITICAL (write2out)
if (bc_timeIncrement(loadcase) < 0.0_pReal) call IO_error(error_ID=34,ext_msg=loadcase_string) ! negative time increment
!$OMP CRITICAL (write2out)
print '(a,f12.6)','Temperature: ',bc_temperature(loadcase) print '(a,f12.6)','Temperature: ',bc_temperature(loadcase)
print '(a,f12.6)','Time: ',bc_timeIncrement(loadcase) print '(a,f12.6)','Time: ',bc_timeIncrement(loadcase)
if (bc_steps(loadcase) < 1_pInt) call IO_error(35,loadcase) ! non-positive increment count !$OMP END CRITICAL (write2out)
if (bc_steps(loadcase) < 1_pInt) call IO_error(error_ID=35,ext_msg=loadcase_string) ! non-positive increment count
!$OMP CRITICAL (write2out)
print '(a,i5)','Increments: ',bc_steps(loadcase) print '(a,i5)','Increments: ',bc_steps(loadcase)
if (bc_frequency(loadcase) < 1_pInt) call IO_error(36,loadcase) ! non-positive result frequency !$OMP END CRITICAL (write2out)
if (bc_frequency(loadcase) < 1_pInt) call IO_error(error_ID=36,ext_msg=loadcase_string) ! non-positive result frequency
!$OMP CRITICAL (write2out)
print '(a,i5)','Freq. of Output: ',bc_frequency(loadcase) print '(a,i5)','Freq. of Output: ',bc_frequency(loadcase)
!$OMP END CRITICAL (write2out)
enddo enddo
!$OMP END CRITICAL (write2out)
ielem = 0_pInt ielem = 0_pInt
c_current = 0.0_pReal c_current = 0.0_pReal
@ -463,7 +480,7 @@ program DAMASK_spectral
#ifdef _OPENMP #ifdef _OPENMP
if(DAMASK_NumThreadsInt>0_pInt) then if(DAMASK_NumThreadsInt>0_pInt) then
call dfftw_init_threads(ierr) call dfftw_init_threads(ierr)
if(ierr == 0_pInt) call IO_error(104,ierr) if(ierr == 0_pInt) call IO_error(error_ID=104)
call dfftw_plan_with_nthreads(DAMASK_NumThreadsInt) call dfftw_plan_with_nthreads(DAMASK_NumThreadsInt)
endif endif
#endif #endif
@ -513,9 +530,9 @@ program DAMASK_spectral
write(538), 'logscale', bc_logscale ! one entry per loadcase (0: linear, 1: log) write(538), 'logscale', bc_logscale ! one entry per loadcase (0: linear, 1: log)
write(538), 'frequencies', bc_frequency ! one entry per loadcase write(538), 'frequencies', bc_frequency ! one entry per loadcase
write(538), 'times', bc_timeIncrement ! one entry per loadcase write(538), 'times', bc_timeIncrement ! one entry per loadcase
bc_timeIncrement(1)= bc_timeIncrement(1) + 1_pInt bc_steps(1)= bc_steps(1) + 1_pInt
write(538), 'increments', bc_timeIncrement ! one entry per loadcase write(538), 'increments', bc_steps ! one entry per loadcase ToDo: rename keyword to steps
bc_timeIncrement(1)= bc_timeIncrement(1) - 1_pInt bc_steps(1)= bc_steps(1) - 1_pInt
write(538), 'startingIncrement', writtenOutCounter write(538), 'startingIncrement', writtenOutCounter
write(538), 'eoh' ! end of header write(538), 'eoh' ! end of header
write(538), materialpoint_results(:,1,:) ! initial (non-deformed) results write(538), materialpoint_results(:,1,:) ! initial (non-deformed) results
@ -613,7 +630,7 @@ program DAMASK_spectral
c_reduced(k,j) = c_prev99(n,m) c_reduced(k,j) = c_prev99(n,m)
endif; enddo; endif; enddo endif; enddo; endif; enddo
call math_invert(size_reduced, c_reduced, s_reduced, i, errmatinv) ! invert reduced stiffness call math_invert(size_reduced, c_reduced, s_reduced, i, errmatinv) ! invert reduced stiffness
if(errmatinv) call IO_error(800) if(errmatinv) call IO_error(error_ID=800)
s_prev99 = 0.0_pReal ! build full compliance s_prev99 = 0.0_pReal ! build full compliance
k = 0_pInt k = 0_pInt
do n = 1,9 do n = 1,9
@ -633,17 +650,21 @@ program DAMASK_spectral
(err_div > err_div_tol .or. & (err_div > err_div_tol .or. &
err_stress > err_stress_tol)) err_stress > err_stress_tol))
iter = iter + 1_pInt iter = iter + 1_pInt
!$OMP CRITICAL (write2out)
print '(A)', '************************************************************' print '(A)', '************************************************************'
print '(3(A,I5.5,tr2)A)', '**** Loadcase = ',loadcase, 'Step = ',step, 'Iteration = ',iter,'****' print '(3(A,I5.5,tr2)A)', '**** Loadcase = ',loadcase, 'Step = ',step, 'Iteration = ',iter,'****'
print '(A)', '************************************************************' print '(A)', '************************************************************'
!$OMP END CRITICAL (write2out)
workfft = 0.0_pReal ! needed because of the padding for FFTW workfft = 0.0_pReal ! needed because of the padding for FFTW
!************************************************************* !*************************************************************
do n = 1,3; do m = 1,3 do n = 1,3; do m = 1,3
defgrad_av(m,n) = sum(defgrad(:,:,:,m,n)) * wgt defgrad_av(m,n) = sum(defgrad(:,:,:,m,n)) * wgt
enddo; enddo enddo; enddo
!$OMP CRITICAL (write2out)
print '(a,/,3(3(f12.7,x)/))', 'Deformation Gradient:',math_transpose3x3(defgrad_av) print '(a,/,3(3(f12.7,x)/))', 'Deformation Gradient:',math_transpose3x3(defgrad_av)
print '(A,/)', '== Update Stress Field (Constitutive Evaluation P(F)) ======' print '(A,/)', '== Update Stress Field (Constitutive Evaluation P(F)) ======'
!$OMP END CRITICAL (write2out)
ielem = 0_pInt ielem = 0_pInt
do k = 1, resolution(3); do j = 1, resolution(2); do i = 1, resolution(1) do k = 1, resolution(3); do j = 1, resolution(2); do i = 1, resolution(1)
ielem = ielem + 1 ielem = ielem + 1
@ -670,6 +691,8 @@ program DAMASK_spectral
do n = 1,3; do m = 1,3 do n = 1,3; do m = 1,3
pstress_av(m,n) = sum(workfft(1:resolution(1),:,:,m,n)) * wgt pstress_av(m,n) = sum(workfft(1:resolution(1),:,:,m,n)) * wgt
enddo; enddo enddo; enddo
!$OMP CRITICAL (write2out)
print '(a,/,3(3(f12.7,x)/))', 'Piola-Kirchhoff Stress / MPa: ',math_transpose3x3(pstress_av)/1.e6 print '(a,/,3(3(f12.7,x)/))', 'Piola-Kirchhoff Stress / MPa: ',math_transpose3x3(pstress_av)/1.e6
err_stress_tol = 0.0_pReal err_stress_tol = 0.0_pReal
@ -686,6 +709,7 @@ program DAMASK_spectral
print '(a,x,f12.7,/)' , 'Determinant of Deformation Aim: ', math_det3x3(defgradAim) print '(a,x,f12.7,/)' , 'Determinant of Deformation Aim: ', math_det3x3(defgradAim)
endif endif
print '(A,/)', '== Calculating Equilibrium Using Spectral Method ===========' print '(A,/)', '== Calculating Equilibrium Using Spectral Method ==========='
!$OMP END CRITICAL (write2out)
call dfftw_execute_dft_r2c(fftw_plan(1),workfft,workfft) ! FFT of pstress call dfftw_execute_dft_r2c(fftw_plan(1),workfft,workfft) ! FFT of pstress
p_hat_avg = sqrt(maxval (math_eigenvalues3x3(math_mul33x33(workfft(1,1,1,1:3,1:3),& ! L_2 norm of average stress in fourier space, p_hat_avg = sqrt(maxval (math_eigenvalues3x3(math_mul33x33(workfft(1,1,1,1:3,1:3),& ! L_2 norm of average stress in fourier space,
@ -744,8 +768,9 @@ program DAMASK_spectral
do m = 1,3; do n = 1,3 do m = 1,3; do n = 1,3
defgrad(:,:,:,m,n) = defgrad(:,:,:,m,n) + (defgradAim_lab(m,n) - defgrad_av(m,n)) ! anticipated target minus current state defgrad(:,:,:,m,n) = defgrad(:,:,:,m,n) + (defgradAim_lab(m,n) - defgrad_av(m,n)) ! anticipated target minus current state
enddo; enddo enddo; enddo
!$OMP CRITICAL (write2out)
print '(2(a,E10.5)/)', 'Error Divergence = ',err_div, ', Tol. = ', err_div_tol print '(2(a,E10.5)/)', 'Error Divergence = ',err_div, ', Tol. = ', err_div_tol
!$OMP END CRITICAL (write2out)
enddo ! end looping when convergency is achieved enddo ! end looping when convergency is achieved
@ -755,19 +780,23 @@ program DAMASK_spectral
write(538), materialpoint_results(:,1,:) ! write result to file write(538), materialpoint_results(:,1,:) ! write result to file
writtenOutCounter = writtenOutCounter + 1_pInt writtenOutCounter = writtenOutCounter + 1_pInt
endif endif
!$OMP CRITICAL (write2out)
if(err_div<=err_div_tol .and. err_stress<=err_stress_tol) then if(err_div<=err_div_tol .and. err_stress<=err_stress_tol) then
print '(2(A,I5.5),A,/)', '== Step = ',step, ' of Loadcase = ',loadcase, ' Converged ==============' print '(2(A,I5.5),A,/)', '== Step = ',step, ' of Loadcase = ',loadcase, ' Converged =============='
else else
print '(2(A,I5.5),A,/)', '== Step = ',step, ' of Loadcase = ',loadcase, ' NOT Converged ==========' print '(2(A,I5.5),A,/)', '== Step = ',step, ' of Loadcase = ',loadcase, ' NOT Converged =========='
notConvergedCounter = notConvergedCounter + 1 notConvergedCounter = notConvergedCounter + 1
endif endif
!$OMP END CRITICAL (write2out)
enddo ! end looping over steps in current loadcase enddo ! end looping over steps in current loadcase
deallocate(c_reduced) deallocate(c_reduced)
deallocate(s_reduced) deallocate(s_reduced)
enddo ! end looping over loadcases enddo ! end looping over loadcases
!$OMP CRITICAL (write2out)
print '(A,/)', '############################################################' print '(A,/)', '############################################################'
print '(a,i5.5,a)', 'A Total of ', notConvergedCounter, ' Steps did not Converge!' print '(a,i5.5,a)', 'A Total of ', notConvergedCounter, ' Steps did not Converge!'
print '(a,i5.5,a)', 'A Total of ', writtenOutCounter, ' Steps are written to File!' print '(a,i5.5,a)', 'A Total of ', writtenOutCounter, ' Steps are written to File!'
!$OMP END CRITICAL (write2out)
close(538) close(538)
call dfftw_destroy_plan(fftw_plan(1)); call dfftw_destroy_plan(fftw_plan(2)) call dfftw_destroy_plan(fftw_plan(1)); call dfftw_destroy_plan(fftw_plan(2))

53
code/IO.f90
View File

@ -1123,17 +1123,17 @@ endfunction
! and terminate the Marc run with exit #9xxx ! and terminate the Marc run with exit #9xxx
! in ABAQUS either time step is reduced or execution terminated ! in ABAQUS either time step is reduced or execution terminated
!******************************************************************** !********************************************************************
subroutine IO_error(ID,e,i,g,ext_msg) subroutine IO_error(error_ID,e,i,g,ext_msg)
use prec, only: pInt use prec, only: pInt
implicit none implicit none
integer(pInt), intent(in) :: ID integer(pInt), intent(in) :: error_ID
integer(pInt), optional, intent(in) :: e,i,g integer(pInt), optional, intent(in) :: e,i,g
character(len=*), optional, intent(in) :: ext_msg character(len=*), optional, intent(in) :: ext_msg
character(len=120) msg character(len=1024) msg
select case (ID) select case (error_ID)
case (30) case (30)
msg = 'could not open spectral loadcase' msg = 'could not open spectral loadcase'
case (31) case (31)
@ -1373,12 +1373,12 @@ endfunction
!$OMP CRITICAL (write2out) !$OMP CRITICAL (write2out)
write(6,*) write(6,*)
write(6,'(a38)') '+------------------------------------+' write(6,'(a38)') '+------------------------------------+'
write(6,'(a38)') '+ error +' write(6,'(a38)') '+ error +'
write(6,'(a17,i3,a18)') '+ ',ID,' +' write(6,'(a17,i3,a18)') '+ ',error_ID,' +'
write(6,'(a38)') '+ +' write(6,'(a38)') '+ +'
write(6,'(a2,a)') '+ ',msg write(6,'(a2,a)') '+ ', trim(msg)
if (present(ext_msg)) write(6,*) '+ ',ext_msg if (present(ext_msg)) write(6,'(a2,a)') '+ ', trim(ext_msg)
if (present(e)) then if (present(e)) then
if (present(i) .and. present(g)) then if (present(i) .and. present(g)) then
write(6,'(a13,i6,a4,i2,a7,i4,a2)') '+ at element ',e,' IP ',i,' grain ',g,' +' write(6,'(a13,i6,a4,i2,a7,i4,a2)') '+ at element ',e,' IP ',i,' grain ',g,' +'
@ -1388,7 +1388,7 @@ endfunction
endif endif
write(6,'(a38)') '+------------------------------------+' write(6,'(a38)') '+------------------------------------+'
call flush(6) call flush(6)
call quit(9000+ID) call quit(9000+error_ID)
!$OMP END CRITICAL (write2out) !$OMP END CRITICAL (write2out)
! ABAQUS returns in some cases ! ABAQUS returns in some cases
@ -1399,28 +1399,28 @@ endfunction
!******************************************************************** !********************************************************************
! write warning statements to standard out ! write warning statements to standard out
!******************************************************************** !********************************************************************
subroutine IO_warning(ID,e,i,g,ext_msg) subroutine IO_warning(warning_ID,e,i,g,ext_msg)
use prec, only: pInt use prec, only: pInt
implicit none implicit none
integer(pInt), intent(in) :: ID integer(pInt), intent(in) :: warning_ID
integer(pInt), optional, intent(in) :: e,i,g integer(pInt), optional, intent(in) :: e,i,g
character(len=*), optional, intent(in) :: ext_msg character(len=*), optional, intent(in) :: ext_msg
character(len=80) msg character(len=1024) msg
select case (ID) select case (warning_ID)
case (33) case (33_pInt)
msg = 'cannot guess along trajectory for first step of first loadcase' msg = 'cannot guess along trajectory for first step of first loadcase'
case (101) case (101_pInt)
msg = '+ crystallite debugging off... +' msg = '+ crystallite debugging off... +'
case (600) case (600_pInt)
msg = '+ crystallite responds elastically +' msg = '+ crystallite responds elastically +'
case (601) case (601_pInt)
msg = '+ stiffness close to zero +' msg = '+ stiffness close to zero +'
case (650) case (650_pInt)
msg = '+ polar decomposition failed +' msg = '+ polar decomposition failed +'
case (700) case (700_pInt)
msg = '+ unknown crystal symmetry +' msg = '+ unknown crystal symmetry +'
case default case default
msg = '+ unknown warning number... +' msg = '+ unknown warning number... +'
@ -1428,11 +1428,12 @@ endfunction
!$OMP CRITICAL (write2out) !$OMP CRITICAL (write2out)
write(6,*) write(6,*)
write(6,'(a38)') '+------------------------------------+' write(6,'(a38)') '+------------------------------------+'
write(6,'(a38)') '+ warning +' write(6,'(a38)') '+ warning +'
write(6,'(a38)') '+ +' write(6,'(a38)') '+ +'
write(6,'(a38)') msg write(6,'(a17,i3,a18)') '+ ',warning_ID,' +'
if (present(ext_msg)) write(6,*) '+ ',ext_msg write(6,'(a2,a)') '+ ', trim(msg)
if (present(ext_msg)) write(6,'(a2,a)') '+ ', trim(ext_msg)
if (present(e)) then if (present(e)) then
if (present(i)) then if (present(i)) then
if (present(g)) then if (present(g)) then

48
code/mesh.f90
View File

@ -333,7 +333,7 @@
parallelExecution = (parallelExecution .and. (mesh_Nelems == mesh_NcpElems)) ! plus potential killer from non-local constitutive parallelExecution = (parallelExecution .and. (mesh_Nelems == mesh_NcpElems)) ! plus potential killer from non-local constitutive
else else
call IO_error(101) ! cannot open input file call IO_error(error_ID=101) ! cannot open input file
endif endif
FEsolving_execElem = (/1,mesh_NcpElems/) FEsolving_execElem = (/1,mesh_NcpElems/)
@ -1468,7 +1468,7 @@ enddo
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,myPos,1) + 1_pInt headerLength = IO_intValue(line,myPos,1) + 1_pInt
else else
call IO_error(42) call IO_error(error_ID=42)
endif endif
rewind(myUnit) rewind(myUnit)
@ -1583,8 +1583,8 @@ enddo
endif endif
enddo enddo
620 if (mesh_Nnodes < 2) call IO_error(900) 620 if (mesh_Nnodes < 2) call IO_error(error_ID=900)
if (mesh_Nelems == 0) call IO_error(901) if (mesh_Nelems == 0) call IO_error(error_ID=901)
endsubroutine endsubroutine
@ -1664,7 +1664,7 @@ enddo
enddo enddo
620 continue 620 continue
if (mesh_NelemSets == 0) call IO_error(902) if (mesh_NelemSets == 0) call IO_error(error_ID=902)
endsubroutine endsubroutine
@ -1706,7 +1706,7 @@ enddo
mesh_Nmaterials = mesh_Nmaterials + 1_pInt mesh_Nmaterials = mesh_Nmaterials + 1_pInt
enddo enddo
620 if (mesh_Nmaterials == 0) call IO_error(903) 620 if (mesh_Nmaterials == 0) call IO_error(error_ID=903)
endsubroutine endsubroutine
@ -1811,7 +1811,7 @@ enddo
endselect endselect
enddo enddo
620 if (mesh_NcpElems == 0) call IO_error(906) 620 if (mesh_NcpElems == 0) call IO_error(error_ID=906)
endsubroutine endsubroutine
@ -1899,7 +1899,7 @@ enddo
640 do i = 1,elemSet 640 do i = 1,elemSet
! write(6,*)'elemSetName: ',mesh_nameElemSet(i) ! write(6,*)'elemSetName: ',mesh_nameElemSet(i)
! write(6,*)'elems in Elset',mesh_mapElemSet(:,i) ! write(6,*)'elems in Elset',mesh_mapElemSet(:,i)
if (mesh_mapElemSet(1,i) == 0) call IO_error(ID=904,ext_msg=mesh_nameElemSet(i)) if (mesh_mapElemSet(1,i) == 0) call IO_error(error_ID=904,ext_msg=mesh_nameElemSet(i))
enddo enddo
endsubroutine endsubroutine
@ -1961,11 +1961,11 @@ enddo
endif endif
enddo enddo
620 if (count==0) call IO_error(905) 620 if (count==0) call IO_error(error_ID=905)
do i=1,count do i=1,count
! write(6,*)'name of materials: ',i,mesh_nameMaterial(i) ! write(6,*)'name of materials: ',i,mesh_nameMaterial(i)
! write(6,*)'name of elemSets: ',i,mesh_mapMaterial(i) ! write(6,*)'name of elemSets: ',i,mesh_mapMaterial(i)
if (mesh_nameMaterial(i)=='' .or. mesh_mapMaterial(i)=='') call IO_error(905) if (mesh_nameMaterial(i)=='' .or. mesh_mapMaterial(i)=='') call IO_error(error_ID=905)
enddo enddo
endsubroutine endsubroutine
@ -2097,7 +2097,7 @@ enddo
650 call qsort(mesh_mapFEtoCPnode,1,size(mesh_mapFEtoCPnode,2)) 650 call qsort(mesh_mapFEtoCPnode,1,size(mesh_mapFEtoCPnode,2))
if (size(mesh_mapFEtoCPnode) == 0) call IO_error(908) if (size(mesh_mapFEtoCPnode) == 0) call IO_error(error_ID=908)
endsubroutine endsubroutine
@ -2227,7 +2227,7 @@ enddo
660 call qsort(mesh_mapFEtoCPelem,1,size(mesh_mapFEtoCPelem,2)) ! should be mesh_NcpElems 660 call qsort(mesh_mapFEtoCPelem,1,size(mesh_mapFEtoCPelem,2)) ! should be mesh_NcpElems
if (size(mesh_mapFEtoCPelem) < 2) call IO_error(907) if (size(mesh_mapFEtoCPelem) < 2) call IO_error(error_ID=907)
endsubroutine endsubroutine
@ -2347,7 +2347,7 @@ subroutine mesh_marc_count_cpSizes (myUnit)
IO_lc(IO_stringValue(line,myPos,2)) /= 'response' ) & IO_lc(IO_stringValue(line,myPos,2)) /= 'response' ) &
) then ) then
t = FE_mapElemtype(IO_extractValue(IO_lc(IO_stringValue(line,myPos,2)),'type')) ! remember elem type t = FE_mapElemtype(IO_extractValue(IO_lc(IO_stringValue(line,myPos,2)),'type')) ! remember elem type
if (t==0) call IO_error(ID=910,ext_msg='mesh_abaqus_count_cpSizes') if (t==0) call IO_error(error_ID=910,ext_msg='mesh_abaqus_count_cpSizes')
count = IO_countDataLines(myUnit) count = IO_countDataLines(myUnit)
do i = 1,count do i = 1,count
backspace(myUnit) backspace(myUnit)
@ -2410,7 +2410,7 @@ subroutine mesh_marc_count_cpSizes (myUnit)
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,myPos,1) + 1_pInt headerLength = IO_intValue(line,myPos,1) + 1_pInt
else else
call IO_error(42) call IO_error(error_ID=42)
endif endif
rewind(myUnit) rewind(myUnit)
@ -2449,9 +2449,9 @@ subroutine mesh_marc_count_cpSizes (myUnit)
! --- sanity checks --- ! --- sanity checks ---
if ((.not. gotDimension) .or. (.not. gotResolution)) call IO_error(42) if ((.not. gotDimension) .or. (.not. gotResolution)) call IO_error(error_ID=42)
if ((a < 1) .or. (b < 1) .or. (c < 0)) call IO_error(43) ! 1_pInt is already added if ((a < 1) .or. (b < 1) .or. (c < 0)) call IO_error(error_ID=43) ! 1_pInt is already added
if ((x <= 0.0_pReal) .or. (y <= 0.0_pReal) .or. (z <= 0.0_pReal)) call IO_error(44) if ((x <= 0.0_pReal) .or. (y <= 0.0_pReal) .or. (z <= 0.0_pReal)) call IO_error(error_ID=44)
forall (n = 0:mesh_Nnodes-1) forall (n = 0:mesh_Nnodes-1)
mesh_node0(1,n+1) = x * dble(mod(n,a) / (a-1.0_pReal)) mesh_node0(1,n+1) = x * dble(mod(n,a) / (a-1.0_pReal))
@ -2561,7 +2561,7 @@ subroutine mesh_marc_count_cpSizes (myUnit)
endif endif
enddo enddo
670 if (size(mesh_node0,2) /= mesh_Nnodes) call IO_error(909) 670 if (size(mesh_node0,2) /= mesh_Nnodes) call IO_error(error_ID=909)
mesh_node = mesh_node0 mesh_node = mesh_node0
endsubroutine endsubroutine
@ -2597,7 +2597,7 @@ subroutine mesh_marc_count_cpSizes (myUnit)
if (keyword(1:4) == 'head') then if (keyword(1:4) == 'head') then
headerLength = IO_intValue(line,myPos,1) + 1_pInt headerLength = IO_intValue(line,myPos,1) + 1_pInt
else else
call IO_error(42) call IO_error(error_ID=42)
endif endif
rewind(myUnit) rewind(myUnit)
@ -2774,7 +2774,7 @@ subroutine mesh_marc_count_cpSizes (myUnit)
IO_lc(IO_stringValue(line,myPos,2)) /= 'response' ) & IO_lc(IO_stringValue(line,myPos,2)) /= 'response' ) &
) then ) then
t = FE_mapElemtype(IO_extractValue(IO_lc(IO_stringValue(line,myPos,2)),'type')) ! remember elem type t = FE_mapElemtype(IO_extractValue(IO_lc(IO_stringValue(line,myPos,2)),'type')) ! remember elem type
if (t==0) call IO_error(ID=910,ext_msg='mesh_abaqus_build_elements') if (t==0) call IO_error(error_ID=910,ext_msg='mesh_abaqus_build_elements')
count = IO_countDataLines(myUnit) count = IO_countDataLines(myUnit)
do i = 1,count do i = 1,count
backspace(myUnit) backspace(myUnit)
@ -3331,13 +3331,13 @@ character(len=64) fmt
integer(pInt) i,e,n,f,t integer(pInt) i,e,n,f,t
if (mesh_maxValStateVar(1) < 1_pInt) call IO_error(110) ! no homogenization specified if (mesh_maxValStateVar(1) < 1_pInt) call IO_error(error_ID=110) ! no homogenization specified
if (mesh_maxValStateVar(2) < 1_pInt) call IO_error(120) ! no microstructure specified if (mesh_maxValStateVar(2) < 1_pInt) call IO_error(error_ID=120) ! no microstructure specified
allocate (mesh_HomogMicro(mesh_maxValStateVar(1),mesh_maxValStateVar(2))); mesh_HomogMicro = 0_pInt allocate (mesh_HomogMicro(mesh_maxValStateVar(1),mesh_maxValStateVar(2))); mesh_HomogMicro = 0_pInt
do e = 1,mesh_NcpElems do e = 1,mesh_NcpElems
if (mesh_element(3,e) < 1_pInt) call IO_error(110,e) ! no homogenization specified if (mesh_element(3,e) < 1_pInt) call IO_error(error_ID=110,e=e) ! no homogenization specified
if (mesh_element(4,e) < 1_pInt) call IO_error(120,e) ! no microstructure specified if (mesh_element(4,e) < 1_pInt) call IO_error(error_ID=120,e=e) ! no microstructure specified
mesh_HomogMicro(mesh_element(3,e),mesh_element(4,e)) = & mesh_HomogMicro(mesh_element(3,e),mesh_element(4,e)) = &
mesh_HomogMicro(mesh_element(3,e),mesh_element(4,e)) + 1 ! count combinations of homogenization and microstructure mesh_HomogMicro(mesh_element(3,e),mesh_element(4,e)) + 1 ! count combinations of homogenization and microstructure
enddo enddo

BIN
documentation/SpectralMethod/Parallel Fourier Transform.pdf Normal file
View File

Binary file not shown.

2
processing/pre/FromEBSD/SpectralMethodFromEBDS
View File

@ -37,6 +37,6 @@ for line in file_in:
dimY = float(re.findall('\S*',str(line))[8]) dimY = float(re.findall('\S*',str(line))[8])
dimZ = min(dimX/resX,dimY/resY) dimZ = min(dimX/resX,dimY/resY)
file_out3 = open(sys.argv[1]+'.geom','w') file_out3 = open(sys.argv[1]+'.geom','w')
file_out3.write('resolution x '+str(resX)+' y '+str(resY)+' z 1 \ndimension a '+str(dimX)+' b '+str(dimY)+' c '+str(dimZ)+'\nhomogenization 1\n') file_out3.write('resolution a '+str(resX)+' b '+str(resY)+' c 1 \ndimension x '+str(dimX)+' y '+str(dimY)+' z '+str(dimZ)+'\nhomogenization 1\n')
for x in xrange(resX*resY): for x in xrange(resX*resY):
file_out3.write(str(x+1)+'\n') file_out3.write(str(x+1)+'\n')