From 649ec1fca8e5033f4d1172b2bf1e46f2476207ac Mon Sep 17 00:00:00 2001 From: Luv Sharma Date: Fri, 19 Sep 2014 17:59:06 +0000 Subject: [PATCH] Finished applying newstate to Homogenization. --- code/CPFEM.f90 | 89 ++++------- code/Makefile | 1 + code/constitutive.f90 | 7 - code/homogenization.f90 | 132 ++-------------- code/homogenization_RGC.f90 | 241 ++---------------------------- code/homogenization_isostrain.f90 | 77 ++-------- code/homogenization_none.f90 | 26 ++-- code/prec.f90 | 16 +- 8 files changed, 87 insertions(+), 502 deletions(-) diff --git a/code/CPFEM.f90 b/code/CPFEM.f90 index b9c34d19c..69a7ec95c 100644 --- a/code/CPFEM.f90 +++ b/code/CPFEM.f90 @@ -147,12 +147,11 @@ subroutine CPFEM_init use material, only: & homogenization_maxNgrains, & material_phase, & -#ifdef NEWSTATE homogState, & mappingHomogenization, & -#endif phase_plasticity, & - plasticState + plasticState, & + material_Nhomogenization use crystallite, only: & crystallite_F0, & crystallite_Fp0, & @@ -160,14 +159,9 @@ subroutine CPFEM_init crystallite_dPdF0, & crystallite_Tstar0_v, & crystallite_localPlasticity - use homogenization, only: & -#ifndef NEWSTATE - homogenization_state0, & -#endif - homogenization_sizeState implicit none - integer(pInt) :: i,j,k,l,m,ph + integer(pInt) :: i,j,k,l,m,ph,homog write(6,'(/,a)') ' <<<+- CPFEM init -+>>>' write(6,'(a)') ' $Id$' @@ -214,29 +208,26 @@ subroutine CPFEM_init call IO_read_realFile(777,'convergedStateConst',modelName) m = 0_pInt - readInstances: do ph = 1_pInt, size(phase_plasticity) + readPlasticityInstances: do ph = 1_pInt, size(phase_plasticity) do k = 1_pInt, plasticState(ph)%sizeState do l = 1, size(plasticState(ph)%state0(1,:)) m = m+1_pInt read(777,rec=m) plasticState(ph)%state0(k,l) enddo; enddo - enddo readInstances + enddo readPlasticityInstances close (777) call IO_read_realFile(777,'convergedStateHomog',modelName) m = 0_pInt - do k = 1,mesh_NcpElems; do j = 1,mesh_maxNips - do l = 1,homogenization_sizeState(j,k) - m = m+1_pInt -#ifdef NEWSTATE - read(777,rec=m) homogState(mappingHomogenization(2,j,k))%state0(l,mappingHomogenization(1,j,k)) -#else - read(777,rec=m) homogenization_state0(j,k)%p(l) -#endif - enddo - - enddo; enddo + readHomogInstances: do homog = 1_pInt, material_Nhomogenization + do k = 1_pInt, homogState(homog)%sizeState + do l = 1, size(homogState(homog)%state0(1,:)) + m = m+1_pInt + read(777,rec=m) homogState(homog)%state0(k,l) + enddo; enddo + enddo readHomogInstances close (777) + #if defined(Marc4DAMASK) || defined(Abaqus) call IO_read_realFile(777,'convergeddcsdE',modelName,size(CPFEM_dcsdE)) read (777,rec=1) CPFEM_dcsdE @@ -314,14 +305,12 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature, dt, elFE, ip) microstructure_elemhomo, & plasticState, & damageState, & -#ifdef NEWSTATE homogState, & - mappingHomogenization, & -#endif thermalState, & mappingConstitutive, & material_phase, & - phase_plasticity + phase_plasticity, & + material_Nhomogenization use crystallite, only: & crystallite_partionedF,& crystallite_F0, & @@ -334,13 +323,7 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature, dt, elFE, ip) crystallite_Tstar0_v, & crystallite_Tstar_v, & crystallite_temperature - use homogenization, only: & - homogenization_sizeState, & -#ifndef NEWSTATE - homogenization_state, & - homogenization_state0, & -#endif - materialpoint_F, & + use homogenization, only: & materialpoint_F, & materialpoint_F0, & materialpoint_P, & materialpoint_dPdF, & @@ -378,7 +361,7 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature, dt, elFE, ip) #endif integer(pInt) elCP, & ! crystal plasticity element number - i, j, k, l, m, n, ph + i, j, k, l, m, n, ph, homog logical updateJaco ! flag indicating if JAcobian has to be updated #if defined(Marc4DAMASK) || defined(Abaqus) @@ -433,19 +416,10 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature, dt, elFE, ip) endif endif - !$OMP PARALLEL DO - do k = 1,mesh_NcpElems - do j = 1,mesh_maxNips - if (homogenization_sizeState(j,k) > 0_pInt) & -#ifdef NEWSTATE - homogState(mappingHomogenization(2,j,k))%state0(:,mappingHomogenization(1,j,k)) = & - homogState(mappingHomogenization(2,j,k))%state(:,mappingHomogenization(1,j,k)) ! internal state of homogenization scheme -#else - homogenization_state0(j,k)%p = homogenization_state(j,k)%p ! internal state of homogenization scheme -#endif - enddo - enddo - !$OMP END PARALLEL DO + do homog = 1_pInt, material_Nhomogenization + homogState(homog)%state0 = homogState(homog)%state + enddo + ! * dump the last converged values of each essential variable to a binary file @@ -479,27 +453,24 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature, dt, elFE, ip) call IO_write_jobRealFile(777,'convergedStateConst') m = 0_pInt - writeInstances: do ph = 1_pInt, size(phase_plasticity) + writePlasticityInstances: do ph = 1_pInt, size(phase_plasticity) do k = 1_pInt, plasticState(ph)%sizeState do l = 1, size(plasticState(ph)%state0(1,:)) m = m+1_pInt write(777,rec=m) plasticState(ph)%state0(k,l) enddo; enddo - enddo writeInstances + enddo writePlasticityInstances close (777) call IO_write_jobRealFile(777,'convergedStateHomog') m = 0_pInt - do k = 1,mesh_NcpElems; do j = 1,mesh_maxNips - do l = 1,homogenization_sizeState(j,k) - m = m+1_pInt -#ifdef NEWSTATE - write(777,rec=m) homogState(mappingHomogenization(2,j,k))%state0(l,mappingHomogenization(1,j,k)) -#else - write(777,rec=m) homogenization_state0(j,k)%p(l) -#endif - enddo - enddo; enddo + writeHomogInstances: do homog = 1_pInt, material_Nhomogenization + do k = 1_pInt, homogState(homog)%sizeState + do l = 1, size(homogState(homog)%state0(1,:)) + m = m+1_pInt + write(777,rec=m) homogState(homog)%state0(k,l) + enddo; enddo + enddo writeHomogInstances close (777) #if defined(Marc4DAMASK) || defined(Abaqus) diff --git a/code/Makefile b/code/Makefile index 0140c7d99..60883125a 100644 --- a/code/Makefile +++ b/code/Makefile @@ -458,6 +458,7 @@ FEM_utilities.o: FEM_utilities.f90 \ FEZoo.o: $(wildcard FEZoo.f90) \ IO.o $(IGNORE) $(PREFIX) $(COMPILERNAME) $(COMPILE) -c ../private/FEM/code/FEZoo.f90 $(SUFFIX) + touch FEZoo.o CPFEM.o: CPFEM.f90\ homogenization.o diff --git a/code/constitutive.f90 b/code/constitutive.f90 index c54cd3abf..21a9a1e10 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -215,17 +215,10 @@ subroutine constitutive_init #endif #ifdef TODO -!-------------------------------------------------------------------------------------------------- -! write out state size file - call IO_write_jobIntFile(777,'sizeStateConst', size(constitutive_sizeState)) - write (777,rec=1) constitutive_sizeState - close(777) - !-------------------------------------------------------------------------------------------------- ! report constitutive_maxSizeState = maxval(constitutive_sizeState) constitutive_maxSizeDotState = maxval(constitutive_sizeDotState) - constitutive_maxSizePostResults = maxval(constitutive_sizePostResults) if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) then write(6,'(a32,1x,7(i8,1x))') 'constitutive_state0: ', shape(constitutive_state0) diff --git a/code/homogenization.f90 b/code/homogenization.f90 index 545a96fba..45d601780 100644 --- a/code/homogenization.f90 +++ b/code/homogenization.f90 @@ -9,17 +9,12 @@ module homogenization use prec, only: & pInt, & - pReal, & - p_vec + pReal !-------------------------------------------------------------------------------------------------- ! General variables for the homogenization at a material point implicit none private -#ifndef NEWSTATE - type(p_vec), dimension(:,:), allocatable, public :: & - homogenization_state0 !< pointer array to homogenization state at start of FE increment -#endif real(pReal), dimension(:,:,:,:), allocatable, public :: & materialpoint_F0, & !< def grad of IP at start of FE increment materialpoint_F, & !< def grad of IP to be reached at end of FE increment @@ -28,20 +23,12 @@ module homogenization materialpoint_dPdF !< tangent of first P--K stress at IP real(pReal), dimension(:,:,:), allocatable, public :: & materialpoint_results !< results array of material point -#ifndef NEWSTATE - type(p_vec), dimension(:,:), allocatable, public, protected :: & - homogenization_state !< pointer array to current homogenization state (end of converged time step) -#endif - integer(pInt), dimension(:,:), allocatable, public, protected :: & - homogenization_sizeState !< size of state array per grain integer(pInt), public, protected :: & materialpoint_sizeResults, & homogenization_maxSizePostResults real(pReal), dimension(:,:), allocatable, public, protected :: & materialpoint_heat - type(p_vec), dimension(:,:), allocatable, private :: & - homogenization_subState0 !< pointer array to homogenization state at start of homogenization increment real(pReal), dimension(:,:,:,:), allocatable, private :: & materialpoint_subF0, & !< def grad of IP at beginning of homogenization increment materialpoint_subF !< def grad of IP to be reached at end of homog inc @@ -233,12 +220,6 @@ subroutine homogenization_init() !-------------------------------------------------------------------------------------------------- ! allocate and initialize global variables -#ifndef NEWSTATE - allocate(homogenization_state0(mesh_maxNips,mesh_NcpElems)) - allocate(homogenization_subState0(mesh_maxNips,mesh_NcpElems)) - allocate(homogenization_state(mesh_maxNips,mesh_NcpElems)) -#endif - allocate(homogenization_sizeState(mesh_maxNips,mesh_NcpElems), source=0_pInt) allocate(homogenization_sizePostResults(mesh_maxNips,mesh_NcpElems), source=0_pInt) allocate(materialpoint_heat(mesh_maxNips,mesh_NcpElems), source=0.0_pReal) allocate(materialpoint_dPdF(3,3,3,3,mesh_maxNips,mesh_NcpElems), source=0.0_pReal) @@ -265,50 +246,13 @@ subroutine homogenization_init() InstancePosition(myInstance) = InstancePosition(myInstance)+1_pInt mapping(e,1:4) = [instancePosition(myinstance),myinstance,e,i] #endif - select case(homogenization_type(mesh_element(3,e))) - case (HOMOGENIZATION_none_ID) -#ifdef NEWSTATE - homogenization_sizePostResults(i,e) = homogState(mappingHomogenization(2,i,e))%sizePostResults -#else - homogenization_sizePostResults(i,e) = 0_pInt -#endif - case (HOMOGENIZATION_ISOSTRAIN_ID) -#ifdef NEWSTATE - homogenization_sizePostResults(i,e) = homogState(mappingHomogenization(2,i,e))%sizePostResults -#else - homogenization_sizePostResults(i,e) = homogenization_isostrain_sizePostResults(myInstance) -#endif - case (HOMOGENIZATION_RGC_ID) - if (homogenization_RGC_sizeState(myInstance) > 0_pInt) then -#ifdef NEWSTATE - homogenization_sizeState(i,e) = homogState(mappingHomogenization(2,i,e))%sizeState -#else - allocate(homogenization_state0(i,e)%p(homogenization_RGC_sizeState(myInstance))) - allocate(homogenization_subState0(i,e)%p(homogenization_RGC_sizeState(myInstance))) - allocate(homogenization_state(i,e)%p(homogenization_RGC_sizeState(myInstance))) - homogenization_state0(i,e)%p = 0.0_pReal - homogenization_sizeState(i,e) = homogenization_RGC_sizeState(myInstance) -#endif - endif -#ifdef NEWSTATE - homogenization_sizePostResults(i,e) = homogState(mappingHomogenization(2,i,e))%sizePostResults -#else - homogenization_sizePostResults(i,e) = homogenization_RGC_sizePostResults(myInstance) -#endif - end select + homogenization_sizePostResults(i,e) = homogState(mappingHomogenization(2,i,e))%sizePostResults enddo IpLooping enddo elementLooping #ifdef HDF call HDF5_mappingHomogenization(mapping) #endif -!-------------------------------------------------------------------------------------------------- -! write state size file out - call IO_write_jobIntFile(777,'sizeStateHomog',size(homogenization_sizeState)) - write (777,rec=1) homogenization_sizeState - close(777) - - homogenization_maxSizeState = maxval(homogenization_sizeState) homogenization_maxSizePostResults = maxval(homogenization_sizePostResults) materialpoint_sizeResults = 1 & ! grain count + 1 + homogenization_maxSizePostResults & ! homogSize & homogResult @@ -323,12 +267,11 @@ subroutine homogenization_init() write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" if (iand(debug_level(debug_homogenization), debug_levelBasic) /= 0_pInt) then -#ifndef NEWSTATE +#ifdef TODO write(6,'(a32,1x,7(i8,1x))') 'homogenization_state0: ', shape(homogenization_state0) write(6,'(a32,1x,7(i8,1x))') 'homogenization_subState0: ', shape(homogenization_subState0) write(6,'(a32,1x,7(i8,1x))') 'homogenization_state: ', shape(homogenization_state) #endif - write(6,'(a32,1x,7(i8,1x))') 'homogenization_sizeState: ', shape(homogenization_sizeState) write(6,'(a32,1x,7(i8,1x),/)') 'homogenization_sizePostResults: ', shape(homogenization_sizePostResults) write(6,'(a32,1x,7(i8,1x))') 'materialpoint_dPdF: ', shape(materialpoint_dPdF) write(6,'(a32,1x,7(i8,1x))') 'materialpoint_F0: ', shape(materialpoint_F0) @@ -344,7 +287,6 @@ subroutine homogenization_init() write(6,'(a32,1x,7(i8,1x))') 'materialpoint_converged: ', shape(materialpoint_converged) write(6,'(a32,1x,7(i8,1x),/)') 'materialpoint_doneAndHappy: ', shape(materialpoint_doneAndHappy) write(6,'(a32,1x,7(i8,1x),/)') 'materialpoint_results: ', shape(materialpoint_results) - write(6,'(a32,1x,7(i8,1x))') 'maxSizeState: ', homogenization_maxSizeState write(6,'(a32,1x,7(i8,1x))') 'maxSizePostResults: ', homogenization_maxSizePostResults endif flush(6) @@ -379,10 +321,8 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) plasticState, & damageState, & thermalState, & -#ifdef NEWSTATE homogState, & - mappingHomogenization, & -#endif + mappingHomogenization, & mappingConstitutive, & homogenization_Ngrains @@ -472,14 +412,10 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) materialpoint_converged(i,e) = .false. ! pretend failed step of twice the required size materialpoint_requested(i,e) = .true. ! everybody requires calculation endforall -#ifdef NEWSTATE - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), homogenization_sizeState(i,e) > 0_pInt) & - homogState(mappingHomogenization(2,i,e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(mappingHomogenization(2,i,e))%State0(:,mappingHomogenization(1,i,e)) ! ...internal homogenization state -#else - forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), homogenization_sizeState(i,e) > 0_pInt) & - homogenization_subState0(i,e)%p = homogenization_state0(i,e)%p ! ...internal homogenization state -#endif + forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), & + homogState(mappingHomogenization(2,i,e))%sizeState > 0_pInt) & + homogState(mappingHomogenization(2,i,e))%subState0(:,mappingHomogenization(1,i,e)) = & + homogState(mappingHomogenization(2,i,e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state enddo NiterationHomog = 0_pInt @@ -526,13 +462,9 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) thermalState(mappingConstitutive(2,g,i,e))%partionedState0(:,mappingConstitutive(1,g,i,e)) = & thermalState(mappingConstitutive(2,g,i,e))%state( :,mappingConstitutive(1,g,i,e)) end forall - if (homogenization_sizeState(i,e) > 0_pInt) & -#ifdef NEWSTATE + if (homogState(mappingHomogenization(2,i,e))%sizeState > 0_pInt) & homogState(mappingHomogenization(2,i,e))%subState0(:,mappingHomogenization(1,i,e)) = & homogState(mappingHomogenization(2,i,e))%state( :,mappingHomogenization(1,i,e)) -#else - homogenization_subState0(i,e)%p = homogenization_state(i,e)%p ! ...internal state of homog scheme -#endif materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) ! ...def grad !$OMP FLUSH(materialpoint_subF0) elseif (materialpoint_requested(i,e)) then steppingNeeded ! already at final time (??) @@ -585,13 +517,9 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) thermalState(mappingConstitutive(2,g,i,e))%state( :,mappingConstitutive(1,g,i,e)) = & thermalState(mappingConstitutive(2,g,i,e))%partionedState0(:,mappingConstitutive(1,g,i,e)) end forall - if (homogenization_sizeState(i,e) > 0_pInt) & -#ifdef NEWSTATE + if (homogState(mappingHomogenization(2,i,e))%sizeState > 0_pInt) & homogState(mappingHomogenization(2,i,e))%state( :,mappingHomogenization(1,i,e)) = & homogState(mappingHomogenization(2,i,e))%subState0(:,mappingHomogenization(1,i,e)) -#else - homogenization_state(i,e)%p = homogenization_subState0(i,e)%p ! ...internal state of homog scheme -#endif endif endif converged @@ -803,21 +731,11 @@ subroutine homogenization_partitionDeformation(ip,el) materialpoint_subF(1:3,1:3,ip,el),& el) case (HOMOGENIZATION_RGC_ID) chosenHomogenization -#ifdef NEWSTATE call homogenization_RGC_partitionDeformation(& crystallite_partionedF(1:3,1:3,1:homogenization_maxNgrains,ip,el), & materialpoint_subF(1:3,1:3,ip,el),& ip, & el) -#else - call homogenization_RGC_partitionDeformation(& - crystallite_partionedF(1:3,1:3,1:homogenization_maxNgrains,ip,el), & - materialpoint_subF(1:3,1:3,ip,el),& - homogenization_state(ip,el), & - ip, & - el) -#endif - end select chosenHomogenization end subroutine homogenization_partitionDeformation @@ -852,7 +770,7 @@ function homogenization_updateState(ip,el) case (HOMOGENIZATION_RGC_ID) chosenHomogenization homogenization_updateState = & -#ifdef NEWSTATE + homogenization_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_maxNgrains,ip,el), & crystallite_partionedF(1:3,1:3,1:homogenization_maxNgrains,ip,el), & crystallite_partionedF0(1:3,1:3,1:homogenization_maxNgrains,ip,el),& @@ -861,18 +779,6 @@ function homogenization_updateState(ip,el) crystallite_dPdF(1:3,1:3,1:3,1:3,1:homogenization_maxNgrains,ip,el), & ip, & el) -#else - homogenization_RGC_updateState(homogenization_state(ip,el), & - homogenization_subState0(ip,el), & - crystallite_P(1:3,1:3,1:homogenization_maxNgrains,ip,el), & - crystallite_partionedF(1:3,1:3,1:homogenization_maxNgrains,ip,el), & - crystallite_partionedF0(1:3,1:3,1:homogenization_maxNgrains,ip,el),& - materialpoint_subF(1:3,1:3,ip,el),& - materialpoint_subdt(ip,el), & - crystallite_dPdF(1:3,1:3,1:3,1:3,1:homogenization_maxNgrains,ip,el), & - ip, & - el) -#endif case default chosenHomogenization homogenization_updateState = .true. end select chosenHomogenization @@ -1027,11 +933,11 @@ function field_getMassDensity(ip,el) select case(field_thermal_type(material_homog(ip,el))) case (FIELD_THERMAL_ADIABATIC_ID) - field_getMassDensity = 0.0_pReal + field_getMassDensity = 0.0_pReal case (FIELD_THERMAL_CONDUCTION_ID) do ipc = 1, homogenization_Ngrains(mesh_element(3,el)) - field_getMassDensity = field_getMassDensity + lattice_massDensity(material_phase(ipc,ip,el)) + field_getMassDensity = field_getMassDensity + lattice_massDensity(material_phase(ipc,ip,el)) enddo end select @@ -1215,6 +1121,7 @@ end function field_getDAMAGE subroutine field_putDAMAGE(ip,el,fieldDamageValue) ! naming scheme use material, only: & fieldDamage, & + material_homog, & mappingHomogenization, & field_damage_type, & FIELD_DAMAGE_NONLOCAL_ID @@ -1338,22 +1245,11 @@ function homogenization_postResults(ip,el) materialpoint_P(1:3,1:3,ip,el), & materialpoint_F(1:3,1:3,ip,el)) case (HOMOGENIZATION_RGC_ID) chosenHomogenization -#ifdef NEWSTATE homogenization_postResults = homogenization_RGC_postResults(& ip, & el, & materialpoint_P(1:3,1:3,ip,el), & materialpoint_F(1:3,1:3,ip,el)) - -#else - homogenization_postResults = homogenization_RGC_postResults(& - homogenization_state(ip,el),& - ip, & - el, & - materialpoint_P(1:3,1:3,ip,el), & - materialpoint_F(1:3,1:3,ip,el)) - -#endif end select chosenHomogenization end function homogenization_postResults diff --git a/code/homogenization_RGC.f90 b/code/homogenization_RGC.f90 index 13a0e0ec3..8d3e25b6a 100644 --- a/code/homogenization_RGC.f90 +++ b/code/homogenization_RGC.f90 @@ -21,10 +21,8 @@ module homogenization_RGC homogenization_RGC_sizePostResult character(len=64), dimension(:,:), allocatable,target, public :: & homogenization_RGC_output ! name of each post result output -#ifdef NEWSTATE integer(pInt), dimension(:), allocatable, private :: & homogenization_RGC_Noutput !< number of outputs per homog instance -#endif integer(pInt), dimension(:,:), allocatable, private :: & homogenization_RGC_Ngrains real(pReal), dimension(:,:), allocatable, private :: & @@ -49,7 +47,7 @@ module homogenization_RGC avgfirstpiola_ID end enum integer(kind(undefined_ID)), dimension(:,:), allocatable, private :: & - homogenization_RGC_outputID !< ID of each post result output + homogenization_RGC_outputID !< ID of each post result output public :: & homogenization_RGC_init, & @@ -103,17 +101,15 @@ subroutine homogenization_RGC_init(fileUnit) use material implicit none - integer(pInt), intent(in) :: fileUnit !< file pointer to material configuration + integer(pInt), intent(in) :: fileUnit !< file pointer to material configuration integer(pInt), parameter :: MAXNCHUNKS = 4_pInt integer(pInt), dimension(1_pInt+2_pInt*MAXNCHUNKS) :: positions -#ifdef NEWSTATE integer :: & homog, & NofMyHomog, & o, & instance, & sizeHState -#endif integer(pInt) :: section=0_pInt, maxNinstance, i,j,e, output=-1_pInt, mySize, myInstance character(len=65536) :: & tag = '', & @@ -126,15 +122,11 @@ subroutine homogenization_RGC_init(fileUnit) maxNinstance = int(count(homogenization_type == HOMOGENIZATION_RGC_ID),pInt) if (maxNinstance == 0_pInt) return -#ifdef NEWSTATE if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0_pInt) & write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance -#endif allocate(homogenization_RGC_sizeState(maxNinstance), source=0_pInt) allocate(homogenization_RGC_sizePostResults(maxNinstance), source=0_pInt) -#ifdef NEWSTATE allocate(homogenization_RGC_Noutput(maxNinstance), source=0_pInt) -#endif allocate(homogenization_RGC_Ngrains(3,maxNinstance), source=0_pInt) allocate(homogenization_RGC_ciAlpha(maxNinstance), source=0.0_pReal) allocate(homogenization_RGC_xiAlpha(maxNinstance), source=0.0_pReal) @@ -177,54 +169,34 @@ subroutine homogenization_RGC_init(fileUnit) homogenization_RGC_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) select case(homogenization_RGC_output(output,i)) case('constitutivework') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = constitutivework_ID case('penaltyenergy') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = penaltyenergy_ID case('volumediscrepancy') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = volumediscrepancy_ID case('averagerelaxrate') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = averagerelaxrate_ID case('maximumrelaxrate') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = maximumrelaxrate_ID case('magnitudemismatch') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = magnitudemismatch_ID case('temperature') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = temperature_ID case('ipcoords') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = ipcoords_ID case('avgdefgrad','avgf') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = avgdefgrad_ID case('avgp','avgfirstpiola','avg1stpiola') -#ifdef NEWSTATE homogenization_RGC_Noutput(i) = homogenization_RGC_Noutput(i) + 1_pInt -#endif homogenization_RGC_outputID(output,i) = avgfirstpiola_ID case default call IO_error(105_pInt,ext_msg=IO_stringValue(line,positions,2_pInt)//& @@ -289,7 +261,6 @@ subroutine homogenization_RGC_init(fileUnit) enddo endif !-------------------------------------------------------------------------------------------------- -#ifdef NEWSTATE initializeInstances: do homog = 1_pInt, material_Nhomogenization myHomog: if (homogenization_type(homog) == HOMOGENIZATION_RGC_ID) then NofMyHomog = count(material_homog == homog) @@ -329,44 +300,14 @@ subroutine homogenization_RGC_init(fileUnit) ! allocate state arrays homogState(homog)%sizeState = sizeHState homogState(homog)%sizePostResults = homogenization_RGC_sizePostResults(instance) - allocate(homogState(homog)%state0 ( sizeHState,NofMyHomog), source=0.0_pReal) - allocate(homogState(homog)%subState0 ( sizeHState,NofMyHomog), source=0.0_pReal) - allocate(homogState(homog)%state ( sizeHState,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%state0 (sizeHState,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%subState0(sizeHState,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%state (sizeHState,NofMyHomog), source=0.0_pReal) endif myHomog enddo initializeInstances -#else - do i = 1_pInt,maxNinstance - do j = 1_pInt,maxval(homogenization_Noutput) - select case(homogenization_RGC_outputID(j,i)) - case(temperature_ID,constitutivework_ID,penaltyenergy_ID,volumediscrepancy_ID, & - averagerelaxrate_ID,maximumrelaxrate_ID) - mySize = 1_pInt - case(ipcoords_ID,magnitudemismatch_ID) - mySize = 3_pInt - case(avgdefgrad_ID,avgfirstpiola_ID) - mySize = 9_pInt - case default - mySize = 0_pInt - end select - - outputFound: if (mySize > 0_pInt) then - homogenization_RGC_sizePostResult(j,i) = mySize - homogenization_RGC_sizePostResults(i) = & - homogenization_RGC_sizePostResults(i) + mySize - endif outputFound - enddo - homogenization_RGC_sizeState(i) & - = 3_pInt*(homogenization_RGC_Ngrains(1,i)-1_pInt)*homogenization_RGC_Ngrains(2,i)*homogenization_RGC_Ngrains(3,i) & - + 3_pInt*homogenization_RGC_Ngrains(1,i)*(homogenization_RGC_Ngrains(2,i)-1_pInt)*homogenization_RGC_Ngrains(3,i) & - + 3_pInt*homogenization_RGC_Ngrains(1,i)*homogenization_RGC_Ngrains(2,i)*(homogenization_RGC_Ngrains(3,i)-1_pInt) & - + 8_pInt ! (1) Average constitutive work, (2-4) Overall mismatch, (5) Average penalty energy, - ! (6) Volume discrepancy, (7) Avg relaxation rate component, (8) Max relaxation rate component - enddo - -#endif end subroutine homogenization_RGC_init @@ -374,14 +315,7 @@ end subroutine homogenization_RGC_init !-------------------------------------------------------------------------------------------------- !> @brief partitions the deformation gradient onto the constituents !-------------------------------------------------------------------------------------------------- -#ifdef NEWSTATE subroutine homogenization_RGC_partitionDeformation(F,avgF,ip,el) -#else -subroutine homogenization_RGC_partitionDeformation(F,avgF,state,ip,el) -#endif - - use prec, only: & - p_vec use debug, only: & debug_level, & debug_homogenization, & @@ -391,10 +325,8 @@ subroutine homogenization_RGC_partitionDeformation(F,avgF,state,ip,el) use material, only: & homogenization_maxNgrains, & homogenization_Ngrains,& -#ifdef NEWSTATE homogState, & - mappingHomogenization, & -#endif + mappingHomogenization, & homogenization_typeInstance use FEsolving, only: & theInc,& @@ -403,9 +335,6 @@ subroutine homogenization_RGC_partitionDeformation(F,avgF,state,ip,el) implicit none real(pReal), dimension (3,3,homogenization_maxNgrains), intent(out) :: F !< partioned F per grain real(pReal), dimension (3,3), intent(in) :: avgF !< averaged F -#ifndef NEWSTATE - type(p_vec), intent(in) :: state -#endif integer(pInt), intent(in) :: & ip, & !< integration point number el !< element number @@ -437,11 +366,9 @@ subroutine homogenization_RGC_partitionDeformation(F,avgF,state,ip,el) iGrain3 = homogenization_RGC_grain1to3(iGrain,homID) do iFace = 1_pInt,nFace intFace = homogenization_RGC_getInterface(iFace,iGrain3) ! identifying 6 interfaces of each grain -#ifdef NEWSTATE + aVect = homogenization_RGC_relaxationVector(intFace,homID, ip, el) ! get the relaxation vectors for each interface from global relaxation vector array -#else - aVect = homogenization_RGC_relaxationVector(intFace,state,homID) ! get the relaxation vectors for each interface from global relaxation vector array -#endif + nVect = homogenization_RGC_interfaceNormal(intFace,ip,el) ! get the normal of each interface forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt) & F(i,j,iGrain) = F(i,j,iGrain) + aVect(i)*nVect(j) ! calculating deformation relaxations due to interface relaxation @@ -470,13 +397,7 @@ end subroutine homogenization_RGC_partitionDeformation !> @brief update the internal state of the homogenization scheme and tell whether "done" and ! "happy" with result !-------------------------------------------------------------------------------------------------- -#ifdef NEWSTATE -function homogenization_RGC_updateState( P,F,F0,avgF,dt,dPdF,ip,el) -#else -function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el) -#endif - use prec, only: & - p_vec +function homogenization_RGC_updateState(P,F,F0,avgF,dt,dPdF,ip,el) use debug, only: & debug_level, & debug_homogenization,& @@ -490,10 +411,8 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el use material, only: & homogenization_maxNgrains, & homogenization_typeInstance, & -#ifdef NEWSTATE homogState, & - mappingHomogenization, & -#endif + mappingHomogenization, & homogenization_Ngrains use numerics, only: & absTol_RGC, & @@ -507,10 +426,7 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el refRelaxRate_RGC implicit none -#ifndef NEWSTATE - type(p_vec), intent(inout) :: state !< current state - type(p_vec), intent(in) :: state0 !< initial state -#endif + real(pReal), dimension (3,3,homogenization_maxNgrains), intent(in) :: & P,& !< array of P F,& !< array of F @@ -554,28 +470,19 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el ! allocate the size of the global relaxation arrays/jacobian matrices depending on the size of the cluster allocate(resid(3_pInt*nIntFaceTot), source=0.0_pReal) allocate(tract(nIntFaceTot,3), source=0.0_pReal) -#ifdef NEWSTATE allocate(relax(3_pInt*nIntFaceTot)); relax= homogState(mappingHomogenization(2,ip,el))% & state(1:3_pInt*nIntFaceTot,mappingHomogenization(1,ip,el)) allocate(drelax(3_pInt*nIntFaceTot)); drelax= homogState(mappingHomogenization(2,ip,el))% & state(1:3_pInt*nIntFaceTot,mappingHomogenization(1,ip,el)) - & homogState(mappingHomogenization(2,ip,el))% & state0(1:3_pInt*nIntFaceTot,mappingHomogenization(1,ip,el)) -#else - allocate(relax(3_pInt*nIntFaceTot)); relax=state%p(1:3_pInt*nIntFaceTot) - allocate(drelax(3_pInt*nIntFaceTot)); drelax=state%p(1:3_pInt*nIntFaceTot)-state0%p(1:3_pInt*nIntFaceTot) -#endif !-------------------------------------------------------------------------------------------------- ! debugging the obtained state if (iand(debug_level(debug_homogenization),debug_levelExtensive) /= 0_pInt) then !$OMP CRITICAL (write2out) write(6,'(1x,a30)')'Obtained state: ' do i = 1_pInt,3_pInt*nIntFaceTot -#ifdef NEWSTATE write(6,'(1x,2(e15.8,1x))')homogState(mappingHomogenization(2,ip,el))%state(i,mappingHomogenization(1,ip,el)) -#else - write(6,'(1x,2(e15.8,1x))')state%p(i) -#endif enddo write(6,*)' ' !$OMP END CRITICAL (write2out) @@ -607,7 +514,7 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el !$OMP END CRITICAL (write2out) endif -!!!------------------------------------------------------------------------------------------------ +!------------------------------------------------------------------------------------------------ ! computing the residual stress from the balance of traction at all (interior) interfaces do iNum = 1_pInt,nIntFaceTot faceID = homogenization_RGC_interface1to4(iNum,homID) ! identifying the interface ID in local coordinate system (4-dimensional index) @@ -689,13 +596,8 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el !-------------------------------------------------------------------------------------------------- ! compute/update the state for postResult, i.e., all energy densities computed by time-integration -#ifdef NEWSTATE constitutiveWork = homogState(mappingHomogenization(2,ip,el))%state(3*nIntFaceTot+1,mappingHomogenization(1,ip,el)) penaltyEnergy = homogState(mappingHomogenization(2,ip,el))%state(3*nIntFaceTot+5,mappingHomogenization(1,ip,el)) -#else - constitutiveWork = state%p(3*nIntFaceTot+1) - penaltyEnergy = state%p(3*nIntFaceTot+5) -#endif do iGrain = 1_pInt,homogenization_Ngrains(mesh_element(3,el)) ! time-integration loop for the calculating the work and energy do i = 1_pInt,3_pInt do j = 1_pInt,3_pInt @@ -704,7 +606,6 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el enddo enddo enddo -#ifdef NEWSTATE homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+1,mappingHomogenization(1,ip,el)) = constitutiveWork ! the bulk mechanical/constitutive work homogState(mappingHomogenization(2,ip,el))% & @@ -723,17 +624,6 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+8,mappingHomogenization(1,ip,el)) = maxval(abs(drelax))/dt ! the maximum rate of relaxation vectors -#else - state%p(3*nIntFaceTot+1) = constitutiveWork ! the bulk mechanical/constitutive work - state%p(3*nIntFaceTot+2) = sum(NN(1,:))/real(nGrain,pReal) ! the overall mismatch of all interface normal to e1-direction - state%p(3*nIntFaceTot+3) = sum(NN(2,:))/real(nGrain,pReal) ! the overall mismatch of all interface normal to e2-direction - state%p(3*nIntFaceTot+4) = sum(NN(3,:))/real(nGrain,pReal) ! the overall mismatch of all interface normal to e3-direction - state%p(3*nIntFaceTot+5) = penaltyEnergy ! the overall penalty energy - state%p(3*nIntFaceTot+6) = volDiscrep ! the overall volume discrepancy - state%p(3*nIntFaceTot+7) = sum(abs(drelax))/dt/real(3_pInt*nIntFaceTot,pReal) ! the average rate of relaxation vectors - state%p(3*nIntFaceTot+8) = maxval(abs(drelax))/dt ! the maximum rate of relaxation vectors - -#endif if (iand(debug_level(debug_homogenization),debug_levelExtensive) /= 0_pInt & .and. debug_e == el .and. debug_i == ip) then !$OMP CRITICAL (write2out) @@ -778,7 +668,7 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el endif -!!!------------------------------------------------------------------------------------------------ +!--------------------------------------------------------------------------------------------------- ! construct the global Jacobian matrix for updating the global relaxation vector array when ! convergence is not yet reached ... @@ -848,13 +738,8 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el do ipert = 1_pInt,3_pInt*nIntFaceTot p_relax = relax p_relax(ipert) = relax(ipert) + pPert_RGC ! perturb the relaxation vector -#ifdef NEWSTATE homogState(mappingHomogenization(2,ip,el))%state(1:3*nIntFaceTot,mappingHomogenization(1,ip,el)) = p_relax call homogenization_RGC_grainDeformation(pF,avgF,ip,el) ! compute the grains deformation from perturbed state -#else - state%p(1:3*nIntFaceTot) = p_relax - call homogenization_RGC_grainDeformation(pF,avgF,state,ip,el) ! compute the grains deformation from perturbed state -#endif call homogenization_RGC_stressPenalty(pR,pNN,avgF,pF,ip,el,homID) ! compute stress penalty due to interface mismatch from perturbed state call homogenization_RGC_volumePenalty(pD,volDiscrep,pF,avgF,ip,el) ! compute stress penalty due to volume discrepancy from perturbed state @@ -969,11 +854,7 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el enddo enddo relax = relax + drelax ! Updateing the state variable for the next iteration -#ifdef NEWSTATE homogState(mappingHomogenization(2,ip,el))%state(1:3*nIntFaceTot,mappingHomogenization(1,ip,el)) = relax -#else - state%p(1:3*nIntFaceTot) = relax -#endif if (any(abs(drelax) > maxdRelax_RGC)) then ! Forcing cutback when the incremental change of relaxation vector becomes too large homogenization_RGC_updateState = [.true.,.false.] !$OMP CRITICAL (write2out) @@ -989,11 +870,7 @@ function homogenization_RGC_updateState( state, state0,P,F,F0,avgF,dt,dPdF,ip,el !$OMP CRITICAL (write2out) write(6,'(1x,a30)')'Returned state: ' do i = 1_pInt,3_pInt*nIntFaceTot -#ifdef NEWSTATE write(6,'(1x,2(e15.8,1x))')homogState(mappingHomogenization(2,ip,el))%state(i,mappingHomogenization(1,ip,el)) -#else - write(6,'(1x,2(e15.8,1x))')state%p(i) -#endif enddo write(6,*)' ' flush(6) @@ -1016,10 +893,8 @@ subroutine homogenization_RGC_averageStressAndItsTangent(avgP,dAvgPdAvgF,P,dPdF, use mesh, only: mesh_element use material, only: & homogenization_maxNgrains, & -#ifdef NEWSTATE homogState, & mappingHomogenization, & -#endif homogenization_Ngrains, & homogenization_typeInstance use math, only: math_Plain3333to99 @@ -1064,22 +939,14 @@ end subroutine homogenization_RGC_averageStressAndItsTangent !-------------------------------------------------------------------------------------------------- !> @brief return array of homogenization results for post file inclusion !-------------------------------------------------------------------------------------------------- -#ifdef NEWSTATE pure function homogenization_RGC_postResults(ip,el,avgP,avgF) -#else -pure function homogenization_RGC_postResults(state,ip,el,avgP,avgF) -#endif - use prec, only: & - p_vec use mesh, only: & mesh_element, & mesh_ipCoordinates use material, only: & homogenization_typeInstance,& -#ifdef NEWSTATE homogState, & - mappingHomogenization, & -#endif + mappingHomogenization, & homogenization_Noutput use crystallite, only: & crystallite_temperature @@ -1091,10 +958,7 @@ pure function homogenization_RGC_postResults(state,ip,el,avgP,avgF) real(pReal), dimension(3,3), intent(in) :: & avgP, & !< average stress at material point avgF !< average deformation gradient at material point -#ifndef NEWSTATE - type(p_vec), intent(in) :: & - state ! my State -#endif + integer(pInt) homID,o,c,nIntFaceTot real(pReal), dimension(homogenization_RGC_sizePostResults(homogenization_typeInstance(mesh_element(3,el)))) :: & homogenization_RGC_postResults @@ -1121,63 +985,33 @@ pure function homogenization_RGC_postResults(state,ip,el,avgP,avgF) homogenization_RGC_postResults(c+1_pInt:c+3_pInt) = mesh_ipCoordinates(1:3,ip,el) ! current ip coordinates c = c + 3_pInt case (constitutivework_ID) -#ifdef NEWSTATE homogenization_RGC_postResults(c+1) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+1,mappingHomogenization(1,ip,el)) -#else - homogenization_RGC_postResults(c+1) = state%p(3*nIntFaceTot+1) -#endif c = c + 1_pInt case (magnitudemismatch_ID) -#ifdef NEWSTATE homogenization_RGC_postResults(c+1) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+2,mappingHomogenization(1,ip,el)) homogenization_RGC_postResults(c+2) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+3,mappingHomogenization(1,ip,el)) homogenization_RGC_postResults(c+3) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+4,mappingHomogenization(1,ip,el)) -#else - homogenization_RGC_postResults(c+1) = state%p(3*nIntFaceTot+2) - homogenization_RGC_postResults(c+2) = state%p(3*nIntFaceTot+3) - homogenization_RGC_postResults(c+3) = state%p(3*nIntFaceTot+4) -#endif c = c + 3_pInt case (penaltyenergy_ID) -#ifdef NEWSTATE homogenization_RGC_postResults(c+1) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+5,mappingHomogenization(1,ip,el)) c = c + 1_pInt -#else - homogenization_RGC_postResults(c+1) = state%p(3*nIntFaceTot+5) - c = c + 1_pInt -#endif case (volumediscrepancy_ID) -#ifdef NEWSTATE homogenization_RGC_postResults(c+1) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+6,mappingHomogenization(1,ip,el)) c = c + 1_pInt -#else - homogenization_RGC_postResults(c+1) = state%p(3*nIntFaceTot+6) - c = c + 1_pInt -#endif case (averagerelaxrate_ID) -#ifdef NEWSTATE homogenization_RGC_postResults(c+1) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+7,mappingHomogenization(1,ip,el)) c = c + 1_pInt -#else - homogenization_RGC_postResults(c+1) = state%p(3*nIntFaceTot+7) - c = c + 1_pInt -#endif case (maximumrelaxrate_ID) -#ifdef NEWSTATE homogenization_RGC_postResults(c+1) = homogState(mappingHomogenization(2,ip,el))% & state(3*nIntFaceTot+8,mappingHomogenization(1,ip,el)) c = c + 1_pInt -#else - homogenization_RGC_postResults(c+1) = state%p(3*nIntFaceTot+8) - c = c + 1_pInt -#endif end select enddo @@ -1203,10 +1037,8 @@ subroutine homogenization_RGC_stressPenalty(rPen,nMis,avgF,fDef,ip,el,homID) math_invert33 use material, only: & homogenization_maxNgrains,& -#ifdef NEWSTATE homogState, & mappingHomogenization, & -#endif homogenization_Ngrains use numerics, only: & xSmoo_RGC @@ -1346,10 +1178,6 @@ subroutine homogenization_RGC_volumePenalty(vPen,vDiscrep,fDef,fAvg,ip,el) math_inv33 use material, only: & homogenization_maxNgrains,& -#ifdef NEWSTATE - homogState, & - mappingHomogenization, & -#endif homogenization_Ngrains use numerics, only: & maxVolDiscr_RGC,& @@ -1409,11 +1237,6 @@ function homogenization_RGC_surfaceCorrection(avgF,ip,el) use math, only: & math_invert33, & math_mul33x33 -#ifdef NEWSTATE - use material, only: & - homogState, & - mappingHomogenization -#endif implicit none real(pReal), dimension(3) :: homogenization_RGC_surfaceCorrection @@ -1483,28 +1306,15 @@ end function homogenization_RGC_equivalentModuli !-------------------------------------------------------------------------------------------------- !> @brief collect relaxation vectors of an interface !-------------------------------------------------------------------------------------------------- -#ifdef NEWSTATE function homogenization_RGC_relaxationVector(intFace,homID, ip, el) -#else -function homogenization_RGC_relaxationVector(intFace,state,homID) -#endif - use prec, only: & - p_vec -#ifdef NEWSTATE use material, only: & homogState, & mappingHomogenization -#endif implicit none -#ifdef NEWSTATE integer(pInt), intent(in) :: ip, el -#endif real(pReal), dimension (3) :: homogenization_RGC_relaxationVector integer(pInt), dimension (4), intent(in) :: intFace !< set of interface ID in 4D array (normal and position) -#ifndef NEWSTATE - type(p_vec), intent(in) :: state !< set of global relaxation vectors -#endif integer(pInt), dimension (3) :: nGDim integer(pInt) :: & iNum, & @@ -1515,12 +1325,8 @@ function homogenization_RGC_relaxationVector(intFace,state,homID) homogenization_RGC_relaxationVector = 0.0_pReal nGDim = homogenization_RGC_Ngrains(1:3,homID) iNum = homogenization_RGC_interface4to1(intFace,homID) ! identify the position of the interface in global state array -#ifdef NEWSTATE if (iNum > 0_pInt) homogenization_RGC_relaxationVector = homogState(mappingHomogenization(2,ip,el))% & state((3*iNum-2):(3*iNum),mappingHomogenization(1,ip,el)) ! get the corresponding entries -#else - if (iNum > 0_pInt) homogenization_RGC_relaxationVector = state%p((3*iNum-2):(3*iNum)) ! get the corresponding entries -#endif end function homogenization_RGC_relaxationVector @@ -1734,30 +1540,17 @@ end function homogenization_RGC_interface1to4 !> @brief calculating the grain deformation gradient (the same with ! homogenization_RGC_partionDeformation, but used only for perturbation scheme) !-------------------------------------------------------------------------------------------------- -#ifdef NEWSTATE subroutine homogenization_RGC_grainDeformation(F, avgF, ip, el) -#else -subroutine homogenization_RGC_grainDeformation(F, avgF, state, ip, el) -#endif - use prec, only: & - p_vec use mesh, only: & mesh_element use material, only: & homogenization_maxNgrains,& homogenization_Ngrains, & -#ifdef NEWSTATE - homogState, & - mappingHomogenization, & -#endif homogenization_typeInstance implicit none real(pReal), dimension (3,3,homogenization_maxNgrains), intent(out) :: F !< partioned F per grain real(pReal), dimension (3,3), intent(in) :: avgF !< -#ifndef NEWSTATE - type(p_vec), intent(in) :: state -#endif integer(pInt), intent(in) :: & el, & !< element number ip !< integration point number @@ -1775,11 +1568,7 @@ subroutine homogenization_RGC_grainDeformation(F, avgF, state, ip, el) iGrain3 = homogenization_RGC_grain1to3(iGrain,homID) do iFace = 1_pInt,nFace intFace = homogenization_RGC_getInterface(iFace,iGrain3) -#ifdef NEWSTATE aVect = homogenization_RGC_relaxationVector(intFace,homID, ip, el) -#else - aVect = homogenization_RGC_relaxationVector(intFace,state,homID) -#endif nVect = homogenization_RGC_interfaceNormal(intFace,ip,el) forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt) & F(i,j,iGrain) = F(i,j,iGrain) + aVect(i)*nVect(j) ! effective relaxations diff --git a/code/homogenization_isostrain.f90 b/code/homogenization_isostrain.f90 index 7342b462f..82d109677 100644 --- a/code/homogenization_isostrain.f90 +++ b/code/homogenization_isostrain.f90 @@ -18,10 +18,8 @@ module homogenization_isostrain character(len=64), dimension(:,:), allocatable, target, public :: & homogenization_isostrain_output !< name of each post result output -#ifdef NEWSTATE - integer(pInt), dimension(:), allocatable, private :: & + integer(pInt), dimension(:), allocatable, private :: & homogenization_isostrain_Noutput !< number of outputs per homog instance -#endif integer(pInt), dimension(:), allocatable, private :: & homogenization_isostrain_Ngrains enum, bind(c) @@ -56,14 +54,11 @@ contains subroutine homogenization_isostrain_init(fileUnit) use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) use prec, only: & - pReal, & - pInt -#ifdef NEWSTATE + pReal use debug, only: & debug_HOMOGENIZATION, & debug_level, & debug_levelBasic -#endif use IO use material @@ -75,14 +70,10 @@ subroutine homogenization_isostrain_init(fileUnit) section = 0_pInt, i, j, output, mySize, o integer :: & maxNinstance, & -#ifdef NEWSTATE homog, & - NofMyHomog, & - instance, & - sizeHState, & -#endif - k -! no pInt (stores a system dependen value from 'count' + instance + integer :: & + NofMyHomog ! no pInt (stores a system dependen value from 'count' character(len=65536) :: & tag = '', & line = '' @@ -94,16 +85,13 @@ subroutine homogenization_isostrain_init(fileUnit) maxNinstance = count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID) if (maxNinstance == 0) return -#ifdef NEWSTATE - if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0_pInt) & + + if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0_pInt) & write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance -#endif allocate(homogenization_isostrain_sizePostResults(maxNinstance), source=0_pInt) allocate(homogenization_isostrain_sizePostResult(maxval(homogenization_Noutput),maxNinstance), & source=0_pInt) -#ifdef NEWSTATE - allocate(homogenization_isostrain_Noutput(maxNinstance), source=0_pInt) -#endif + allocate(homogenization_isostrain_Noutput(maxNinstance), source=0_pInt) allocate(homogenization_isostrain_Ngrains(maxNinstance), source=0_pInt) allocate(homogenization_isostrain_mapping(maxNinstance), source=average_ID) allocate(homogenization_isostrain_output(maxval(homogenization_Noutput),maxNinstance)) @@ -116,7 +104,7 @@ subroutine homogenization_isostrain_init(fileUnit) line = IO_read(fileUnit) enddo - parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of homogenization part + parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of homogenization part line = IO_read(fileUnit) if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') then ! stop at next part @@ -143,29 +131,19 @@ subroutine homogenization_isostrain_init(fileUnit) homogenization_isostrain_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) select case(homogenization_isostrain_output(output,i)) case('nconstituents','ngrains') -#ifdef NEWSTATE homogenization_isostrain_Noutput(i) = homogenization_isostrain_Noutput(i) + 1_pInt -#endif homogenization_isostrain_outputID(output,i) = nconstituents_ID case('temperature') -#ifdef NEWSTATE homogenization_isostrain_Noutput(i) = homogenization_isostrain_Noutput(i) + 1_pInt -#endif homogenization_isostrain_outputID(output,i) = temperature_ID case('ipcoords') -#ifdef NEWSTATE homogenization_isostrain_Noutput(i) = homogenization_isostrain_Noutput(i) + 1_pInt -#endif homogenization_isostrain_outputID(output,i) = ipcoords_ID case('avgdefgrad','avgf') -#ifdef NEWSTATE homogenization_isostrain_Noutput(i) = homogenization_isostrain_Noutput(i) + 1_pInt -#endif homogenization_isostrain_outputID(output,i) = avgdefgrad_ID case('avgp','avgfirstpiola','avg1stpiola') -#ifdef NEWSTATE homogenization_isostrain_Noutput(i) = homogenization_isostrain_Noutput(i) + 1_pInt -#endif homogenization_isostrain_outputID(output,i) = avgfirstpiola_ID case default call IO_error(105_pInt,ext_msg=IO_stringValue(line,positions,2_pInt)//& @@ -189,7 +167,6 @@ subroutine homogenization_isostrain_init(fileUnit) endif enddo parsingFile -#ifdef NEWSTATE initializeInstances: do homog = 1_pInt, material_Nhomogenization myHomog: if (homogenization_type(homog) == HOMOGENIZATION_ISOSTRAIN_ID) then NofMyHomog = count(material_homog == homog) @@ -216,43 +193,15 @@ subroutine homogenization_isostrain_init(fileUnit) enddo outputsLoop ! allocate state arrays - sizeHState = 0_pInt - homogState(homog)%sizeState = sizeHState + homogState(homog)%sizeState = 0_pInt homogState(homog)%sizePostResults = homogenization_isostrain_sizePostResults(instance) - allocate(homogState(homog)%state0 ( sizeHState,NofMyHomog), source=0.0_pReal) - allocate(homogState(homog)%subState0 ( sizeHState,NofMyHomog), source=0.0_pReal) - allocate(homogState(homog)%state ( sizeHState,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%state0 (0_pInt,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%subState0(0_pInt,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%state (0_pInt,NofMyHomog), source=0.0_pReal) endif myHomog enddo initializeInstances -#else - - do i = 1,maxNinstance - - do j = 1_pInt,maxval(homogenization_Noutput) - select case(homogenization_isostrain_outputID(j,i)) - case(nconstituents_ID, temperature_ID) - mySize = 1_pInt - case(ipcoords_ID) - mySize = 3_pInt - case(avgdefgrad_ID, avgfirstpiola_ID) - mySize = 9_pInt - case default - mySize = 0_pInt - end select - - outputFound: if (mySize > 0_pInt) then - homogenization_isostrain_sizePostResult(j,i) = mySize - homogenization_isostrain_sizePostResults(i) = & - homogenization_isostrain_sizePostResults(i) + mySize - endif outputFound - enddo - enddo - -#endif - - end subroutine homogenization_isostrain_init diff --git a/code/homogenization_none.f90 b/code/homogenization_none.f90 index c58337719..1c07a9bdc 100644 --- a/code/homogenization_none.f90 +++ b/code/homogenization_none.f90 @@ -29,36 +29,28 @@ subroutine homogenization_none_init() use material implicit none -#ifdef NEWSTATE - integer :: & + integer(pInt) :: & homog, & - NofMyHomog, & - instance, & - sizeHState -#endif + NofMyHomog + write(6,'(/,a)') ' <<<+- homogenization_'//HOMOGENIZATION_NONE_label//' init -+>>>' write(6,'(a)') ' $Id$' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() #include "compilation_info.f90" -#ifdef NEWSTATE initializeInstances: do homog = 1_pInt, material_Nhomogenization myhomog: if (homogenization_type(homog) == HOMOGENIZATION_none_ID) then - NofMyHomog = count(material_homog == homog) -! instance = phase_plasticityInstance(phase) - -! allocate homogenization state arrays - sizeHState = 0_pInt - homogState(homog)%sizeState = sizeHState + NofMyHomog = count(material_homog == homog) + homogState(homog)%sizeState = 0_pInt homogState(homog)%sizePostResults = 0_pInt - allocate(homogState(homog)%state0 ( sizeHState,NofMyHomog), source=0.0_pReal) - allocate(homogState(homog)%subState0 ( sizeHState,NofMyHomog), source=0.0_pReal) - allocate(homogState(homog)%state ( sizeHState,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%state0 (0_pInt,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%subState0(0_pInt,NofMyHomog), source=0.0_pReal) + allocate(homogState(homog)%state (0_pInt,NofMyHomog), source=0.0_pReal) endif myhomog enddo initializeInstances -#endif + end subroutine homogenization_none_init diff --git a/code/prec.f90 b/code/prec.f90 index adbc56d04..03920c67d 100644 --- a/code/prec.f90 +++ b/code/prec.f90 @@ -77,18 +77,12 @@ module prec RK4dotState real(pReal), allocatable, dimension(:,:,:) :: RKCK45dotState end type -#ifdef NEWSTATE - type, public :: hState - integer(pInt) :: sizeState = 0_pInt , & + +#ifdef NEWSTATE + type, public :: tFieldData + integer(pInt) :: sizeField = 0_pInt , & sizePostResults = 0_pInt - real(pReal), allocatable, dimension(:,:) :: state, & ! material points, state size - state0, & - subState0 - end type - type, public :: fState - integer(pInt) :: sizeState = 0_pInt , & - sizePostResults = 0_pInt - real(pReal), allocatable, dimension(:,:) :: state ! material points, state size + real(pReal), allocatable, dimension(:,:) :: field ! material points, state size end type #endif