diff --git a/code/CPFEM.f90 b/code/CPFEM.f90 index ec963c3c5..0e7f78f4b 100644 --- a/code/CPFEM.f90 +++ b/code/CPFEM.f90 @@ -116,7 +116,7 @@ subroutine CPFEM_init() use FEsolving, only: parallelExecution, & symmetricSolver, & restartRead, & - restartJob + FEmodelGeometry use mesh, only: mesh_NcpElems, & mesh_maxNips use material, only: homogenization_maxNgrains, & @@ -149,31 +149,31 @@ subroutine CPFEM_init() write(6,'(a)') '<< CPFEM >> Restored state variables of last converged step from binary files' !$OMP END CRITICAL (write2out) endif - if (IO_read_jobBinaryFile(777,'recordedPhase',restartJob,size(material_phase))) then + if (IO_read_jobBinaryFile(777,'recordedPhase',FEmodelGeometry,size(material_phase))) then read (777,rec=1) material_phase close (777) endif - if (IO_read_jobBinaryFile(777,'convergedF',restartJob,size(crystallite_F0))) then + if (IO_read_jobBinaryFile(777,'convergedF',FEmodelGeometry,size(crystallite_F0))) then read (777,rec=1) crystallite_F0 close (777) endif - if (IO_read_jobBinaryFile(777,'convergedFp',restartJob,size(crystallite_Fp0))) then + if (IO_read_jobBinaryFile(777,'convergedFp',FEmodelGeometry,size(crystallite_Fp0))) then read (777,rec=1) crystallite_Fp0 close (777) endif - if (IO_read_jobBinaryFile(777,'convergedLp',restartJob,size(crystallite_Lp0))) then + if (IO_read_jobBinaryFile(777,'convergedLp',FEmodelGeometry,size(crystallite_Lp0))) then read (777,rec=1) crystallite_Lp0 close (777) endif - if (IO_read_jobBinaryFile(777,'convergeddPdF',restartJob,size(crystallite_dPdF0))) then + if (IO_read_jobBinaryFile(777,'convergeddPdF',FEmodelGeometry,size(crystallite_dPdF0))) then read (777,rec=1) crystallite_dPdF0 close (777) endif - if (IO_read_jobBinaryFile(777,'convergedTstar',restartJob,size(crystallite_Tstar0_v))) then + if (IO_read_jobBinaryFile(777,'convergedTstar',FEmodelGeometry,size(crystallite_Tstar0_v))) then read (777,rec=1) crystallite_Tstar0_v close (777) endif - if (IO_read_jobBinaryFile(777,'convergedStateConst',restartJob)) then + if (IO_read_jobBinaryFile(777,'convergedStateConst',FEmodelGeometry)) then m = 0_pInt do i = 1,homogenization_maxNgrains; do j = 1,mesh_maxNips; do k = 1,mesh_NcpElems do l = 1,size(constitutive_state0(i,j,k)%p) @@ -183,7 +183,7 @@ subroutine CPFEM_init() enddo; enddo; enddo close (777) endif - if (IO_read_jobBinaryFile(777,'convergedStateHomog',restartJob)) then + if (IO_read_jobBinaryFile(777,'convergedStateHomog',FEmodelGeometry)) then m = 0_pInt do k = 1,mesh_NcpElems; do j = 1,mesh_maxNips do l = 1,homogenization_sizeState(j,k) @@ -193,7 +193,7 @@ subroutine CPFEM_init() enddo; enddo close (777) endif - if (IO_read_jobBinaryFile(777,'convergeddcsdE',restartJob,size(CPFEM_dcsdE))) then + if (IO_read_jobBinaryFile(777,'convergeddcsdE',FEmodelGeometry,size(CPFEM_dcsdE))) then read (777,rec=1) CPFEM_dcsdE close (777) endif diff --git a/code/FEsolving.f90 b/code/FEsolving.f90 index 6fe7e4b86..99957e147 100644 --- a/code/FEsolving.f90 +++ b/code/FEsolving.f90 @@ -36,7 +36,7 @@ logical, dimension(:,:), allocatable :: calcMode integer(pInt), dimension(:,:), allocatable :: FEsolving_execIP integer(pInt), dimension(2) :: FEsolving_execElem - character(len=1024) restartJob + character(len=1024) FEmodelGeometry CONTAINS @@ -48,17 +48,20 @@ use prec, only: pInt use debug, only: debug_verbosity + use DAMASK_interface, only: getModelName, FEsolver use IO implicit none integer(pInt), parameter :: fileunit = 222 integer(pInt), parameter :: maxNchunks = 6 + integer(pInt) i integer(pInt), dimension(1+2*maxNchunks) :: positions character(len=64) tag character(len=1024) line + FEmodelGeometry = getModelName() - if (IO_open_inputFile(fileunit)) then + if (IO_open_inputFile(fileunit,FEmodelGeometry)) then rewind(fileunit) do @@ -75,6 +78,11 @@ positions = IO_stringPos(line,maxNchunks) restartWrite = iand(IO_intValue(line,positions,1),1_pInt) > 0_pInt restartRead = iand(IO_intValue(line,positions,1),2_pInt) > 0_pInt + case ('*restart') + do i=2,positions(1) + restartWrite = IO_lc(IO_StringValue(line,positions,i)) == 'write' +! restartRead = IO_lc(IO_StringValue(line,positions,i)) == 'read' + enddo end select enddo else @@ -83,17 +91,29 @@ 100 close(fileunit) - if (restartRead .and. IO_open_logFile(fileunit)) then - rewind(fileunit) - do - read (fileunit,'(a1024)',END=200) line - positions = IO_stringPos(line,maxNchunks) - if ( IO_lc(IO_stringValue(line,positions,1)) == 'restart' .and. & - IO_lc(IO_stringValue(line,positions,2)) == 'file' .and. & - IO_lc(IO_stringValue(line,positions,3)) == 'job' .and. & - IO_lc(IO_stringValue(line,positions,4)) == 'id' ) & - restartJob = IO_StringValue(line,positions,6) - enddo + if (restartRead) then + if(FEsolver == 'Marc' .and. IO_open_logFile(fileunit)) then + rewind(fileunit) + do + read (fileunit,'(a1024)',END=200) line + positions = IO_stringPos(line,maxNchunks) + if ( IO_lc(IO_stringValue(line,positions,1)) == 'restart' .and. & + IO_lc(IO_stringValue(line,positions,2)) == 'file' .and. & + IO_lc(IO_stringValue(line,positions,3)) == 'job' .and. & + IO_lc(IO_stringValue(line,positions,4)) == 'id' ) & + FEmodelGeometry = IO_StringValue(line,positions,6) + enddo + elseif (FEsolver == 'Abaqus' .and. IO_open_jobFile(fileunit, 'com')) then + rewind(fileunit) + do + read (fileunit,'(a1024)',END=200) line + positions = IO_stringPos(line,maxNchunks) +! if ( IO_lc(IO_stringValue(line,positions,?)) == 'oldjob?') & +! FEmodelGeometry = IO_StringValue(line,positions,?) + enddo + else + call IO_error(106) ! cannot open file for old job info + endif endif 200 close(fileunit) @@ -106,7 +126,7 @@ if (debug_verbosity > 0) then write(6,*) 'restart writing: ', restartWrite write(6,*) 'restart reading: ', restartRead - if (restartRead) write(6,*) 'restart Job: ', trim(restartJob) + if (restartRead) write(6,*) 'restart Job: ', trim(FEmodelGeometry) write(6,*) endif !$OMP END CRITICAL (write2out) diff --git a/code/IO.f90 b/code/IO.f90 index ca31ce101..55501a112 100644 --- a/code/IO.f90 +++ b/code/IO.f90 @@ -26,7 +26,8 @@ !--------------------------- ! function IO_abaqus_assembleInputFile ! function IO_open_file(unit,relPath) -! function IO_open_inputFile(unit) +! function IO_open_inputFile(unit, model) +! function IO_open_logFile(unit) ! function IO_hybridIA(Nast,ODFfileName) ! private function hybridIA_reps(dV_V,steps,C) ! function IO_stringPos(line,maxN) @@ -175,27 +176,28 @@ end function ! : changed the function to open *.inp_assembly, which is basically ! the input file without comment lines and possibly assembled includes !******************************************************************** - logical function IO_open_inputFile(unit) + logical function IO_open_inputFile(unit,model) use prec, only: pReal, pInt use DAMASK_interface implicit none integer(pInt), intent(in) :: unit + character(len=*) model IO_open_inputFile = .false. if (FEsolver == 'Abaqus') then open(unit+1,status='old',err=100,& file=trim(getSolverWorkingDirectoryName())//& - trim(getModelName())//InputFileExtension) + trim(model)//InputFileExtension) open(unit,err=100,file=trim(getSolverWorkingDirectoryName())//& - trim(getModelName())//InputFileExtension//'_assembly') + trim(model)//InputFileExtension//'_assembly') IO_open_inputFile = IO_abaqus_assembleInputFile(unit,unit+1) ! strip comments and concatenate any "include"s close(unit+1) else open(unit,status='old',err=100,file=trim(getSolverWorkingDirectoryName())//& - trim(getModelName())//InputFileExtension) + trim(model)//InputFileExtension) IO_open_inputFile = .true. endif @@ -1183,6 +1185,8 @@ endfunction msg = 'FFTW init error' case (105) msg = 'Error reading from ODF file' + case (106) + msg = 'Error reading info on old job' case (110) msg = 'No homogenization specified via State Variable 2' case (120) diff --git a/code/config/material.config b/code/config/material.config index 473289aaa..55a5fe9b7 100644 --- a/code/config/material.config +++ b/code/config/material.config @@ -393,3 +393,6 @@ symmetry orthotropic # or monoclinic [fiber example] # fiber axis in spherical coordinates: alpha crystal system, beta sample system (fiber) alpha1 123 alpha2 123 beta1 12 beta2 45 scatter 15 fraction 0.333 + +[random single crystals] +(random) scatter 0.000 fraction 1.000 \ No newline at end of file diff --git a/code/material.f90 b/code/material.f90 index 4e9838782..81c0acff6 100644 --- a/code/material.f90 +++ b/code/material.f90 @@ -445,7 +445,7 @@ subroutine material_parseTexture(file,myPart) use prec, only: pInt, pReal use IO - use math, only: inRad + use math, only: inRad, math_sampleRandomOri implicit none character(len=*), intent(in) :: myPart @@ -467,8 +467,9 @@ subroutine material_parseTexture(file,myPart) allocate(texture_Ngauss(Nsections)); texture_Ngauss = 0_pInt allocate(texture_Nfiber(Nsections)); texture_Nfiber = 0_pInt - texture_Ngauss = IO_countTagInPart(file,myPart,'(gauss)',Nsections) - texture_Nfiber = IO_countTagInPart(file,myPart,'(fiber)',Nsections) + texture_Ngauss = IO_countTagInPart(file,myPart,'(gauss)', Nsections) + & + IO_countTagInPart(file,myPart,'(random)',Nsections) + texture_Nfiber = IO_countTagInPart(file,myPart,'(fiber)', Nsections) texture_maxNgauss = maxval(texture_Ngauss) texture_maxNfiber = maxval(texture_Nfiber) allocate(texture_Gauss (5,texture_maxNgauss,Nsections)); texture_Gauss = 0.0_pReal @@ -511,6 +512,19 @@ subroutine material_parseTexture(file,myPart) texture_symmetry(section) = 1 end select + case ('(random)') + gauss = gauss + 1 + texture_Gauss(1:3,gauss,section) = math_sampleRandomOri() + do i = 2,4,2 + tag = IO_lc(IO_stringValue(line,positions,i)) + select case (tag) + case('scatter') + texture_Gauss(4,gauss,section) = IO_floatValue(line,positions,i+1)*inRad + case('fraction') + texture_Gauss(5,gauss,section) = IO_floatValue(line,positions,i+1) + end select + enddo + case ('(gauss)') gauss = gauss + 1 do i = 2,10,2 @@ -528,6 +542,7 @@ subroutine material_parseTexture(file,myPart) texture_Gauss(5,gauss,section) = IO_floatValue(line,positions,i+1) end select enddo + case ('(fiber)') fiber = fiber + 1 do i = 2,12,2 diff --git a/code/math.f90 b/code/math.f90 index 96ef2c299..1d358d068 100644 --- a/code/math.f90 +++ b/code/math.f90 @@ -166,6 +166,7 @@ real(pReal), dimension(4,36), parameter :: math_symOperations = & call random_seed(get=randInit) if (debug_verbosity > 0) then !$OMP CRITICAL (write2out) + ! this critical block did cause trouble at IWM write(6,*) 'random seed: ',randInit(1) write(6,*) !$OMP END CRITICAL (write2out) diff --git a/code/mesh.f90 b/code/mesh.f90 index ff8aa83ac..c9274cec3 100644 --- a/code/mesh.f90 +++ b/code/mesh.f90 @@ -265,7 +265,7 @@ use DAMASK_interface use prec, only: pInt use IO, only: IO_error,IO_open_InputFile,IO_abaqus_hasNoPart - use FEsolving, only: parallelExecution, FEsolving_execElem, FEsolving_execIP, calcMode, lastMode + use FEsolving, only: parallelExecution, FEsolving_execElem, FEsolving_execIP, calcMode, lastMode, FEmodelGeometry implicit none @@ -279,9 +279,9 @@ write(6,*) !$OMP END CRITICAL (write2out) - call mesh_build_FEdata() ! --- get properties of the different types of elements + call mesh_build_FEdata() ! --- get properties of the different types of elements - if (IO_open_inputFile(fileUnit)) then ! --- parse info from input file... + if (IO_open_inputFile(fileUnit,FEmodelGeometry)) then ! --- parse info from input file... select case (FEsolver) case ('Spectral')