From 986926aaf2e9647e861fe6f765879b769adf41d7 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 8 May 2014 14:55:19 +0000 Subject: [PATCH] some changes related to new state fixed assigning of wrong output in J2 model --- code/constitutive.f90 | 27 ++++++++++---- code/constitutive_j2.f90 | 58 ++++++++++++++++++++++------- code/constitutive_phenopowerlaw.f90 | 4 +- code/material.f90 | 14 ++++--- code/prec.f90 | 16 +++++++- 5 files changed, 88 insertions(+), 31 deletions(-) diff --git a/code/constitutive.f90 b/code/constitutive.f90 index 28c053e6f..83c4929fd 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -86,7 +86,6 @@ subroutine constitutive_init FE_Nips, & FE_geomtype use material, only: & - material_phase, & material_phase, & material_Nphase, & material_localFileExt, & @@ -111,6 +110,9 @@ subroutine constitutive_init PLASTICITY_PHENOPOWERLAW_label, & PLASTICITY_DISLOTWIN_label, & PLASTICITY_TITANMOD_label, & +#ifdef NEWSTATE + plasticState, & +#endif PLASTICITY_NONLOCAL_label use constitutive_none use constitutive_j2 @@ -130,13 +132,14 @@ subroutine constitutive_init eMax, & !< maximum number of elements phase, & s, & - instance,& + instance,& myNgrains + integer(pInt), dimension(:,:), pointer :: thisSize character(len=64), dimension(:,:), pointer :: thisOutput character(len=32) :: outputName !< name of output, intermediate fix until HDF5 output is ready logical :: knownPlasticity, nonlocalConstitutionPresent -#ifdef HDF +#if defined(HDF) || defined(NEWSTATE) integer(pInt), dimension(:,:,:), allocatable :: mappingConstitutive integer(pInt), dimension(:,:,:), allocatable :: mappingCrystallite integer(pInt), dimension(:), allocatable :: ConstitutivePosition @@ -148,7 +151,6 @@ subroutine constitutive_init #endif nonlocalConstitutionPresent = .false. - !-------------------------------------------------------------------------------------------------- ! parse plasticities from config file if (.not. IO_open_jobFile_stat(FILEUNIT,material_localFileExt)) & ! no local material configuration present... @@ -215,7 +217,10 @@ subroutine constitutive_init cMax = homogenization_maxNgrains iMax = mesh_maxNips eMax = mesh_NcpElems + + +! lumped into new state allocate(constitutive_state0(cMax,iMax,eMax)) allocate(constitutive_partionedState0(cMax,iMax,eMax)) allocate(constitutive_subState0(cMax,iMax,eMax)) @@ -225,6 +230,7 @@ subroutine constitutive_init allocate(constitutive_deltaState(cMax,iMax,eMax)) allocate(constitutive_dotState_backup(cMax,iMax,eMax)) allocate(constitutive_aTolState(cMax,iMax,eMax)) +! not needed anymore for new state allocate(constitutive_sizeDotState(cMax,iMax,eMax), source=0_pInt) allocate(constitutive_sizeState(cMax,iMax,eMax), source=0_pInt) allocate(constitutive_sizePostResults(cMax,iMax,eMax), source=0_pInt) @@ -248,7 +254,7 @@ subroutine constitutive_init end select phase = material_phase(g,i,e) instance = phase_plasticityInstance(phase) -#ifdef HDF +#if defined(HDF) || defined(NEWSTATE) ConstitutivePosition(phase) = ConstitutivePosition(phase)+1_pInt mappingConstitutive(g,e,1:2) = [ConstitutivePosition(phase),phase] #endif @@ -280,7 +286,7 @@ subroutine constitutive_init constitutive_sizeState(g,i,e) = 0_pInt constitutive_sizeDotState(g,i,e) = 0_pInt constitutive_sizePostResults(g,i,e) = 0_pInt - case (PLASTICITY_J2_ID) + case (PLASTICITY_J2_ID) allocate(constitutive_state0(g,i,e)%p(constitutive_j2_sizeState(instance))) allocate(constitutive_partionedState0(g,i,e)%p(constitutive_j2_sizeState(instance))) allocate(constitutive_subState0(g,i,e)%p(constitutive_j2_sizeState(instance))) @@ -439,6 +445,8 @@ subroutine constitutive_init end select enddo #endif + + !-------------------------------------------------------------------------------------------------- ! write out state size file call IO_write_jobIntFile(777,'sizeStateConst', size(constitutive_sizeState)) @@ -537,11 +545,14 @@ subroutine constitutive_microstructure(temperature, Fe, Fp, ipc, ip, el) real(pReal), intent(in), dimension(3,3) :: & Fe, & !< elastic deformation gradient Fp !< plastic deformation gradient - + ! offset = mappingConstitutive(ipc,el,1) ???? + +! OLD constitutive_state(ipc,ip,el) +! NEW plasticState(phase=material_phase(ipc,ip,el))%state(1:myStateSize,offset) select case (phase_plasticity(material_phase(ipc,ip,el))) case (PLASTICITY_DISLOTWIN_ID) - call constitutive_dislotwin_microstructure(temperature,constitutive_state(ipc,ip,el), & + call constitutive_dislotwin_microstructure(temperature,constitutive_state(ipc,ip,el), & ipc,ip,el) case (PLASTICITY_TITANMOD_ID) call constitutive_titanmod_microstructure(temperature,constitutive_state(ipc,ip,el), & diff --git a/code/constitutive_j2.f90 b/code/constitutive_j2.f90 index e0c2b8d70..b1dcae24d 100644 --- a/code/constitutive_j2.f90 +++ b/code/constitutive_j2.f90 @@ -92,10 +92,12 @@ subroutine constitutive_j2_init(fileUnit) #ifdef HDF use hdf5 #endif -use debug, only: & + use debug, only: & debug_level, & debug_constitutive, & debug_levelBasic + use numerics, only: & + numerics_integrator use math, only: & math_Mandel3333to66, & math_Voigt66to3333 @@ -122,8 +124,12 @@ use debug, only: & phase_Noutput, & PLASTICITY_J2_label, & PLASTICITY_J2_ID, & - material_phase,& + material_phase, & +#ifdef NEWSTATE + plasticState, & +#endif MATERIAL_partPhase + use lattice implicit none @@ -136,12 +142,15 @@ use debug, only: & character(len=65536) :: & tag = '', & line = '' + integer(pInt) :: NofMyPhase #ifdef HDF character(len=5) :: & str1 integer(HID_T) :: ID,ID2,ID4 #endif - +#ifdef NEWSTATE + integer(pInt) :: sizeDotState,sizeState +#endif write(6,'(/,a)') ' <<<+- constitutive_'//PLASTICITY_J2_label//' init -+>>>' write(6,'(a)') ' $Id$' write(6,'(a15,a)') ' Current time: ',IO_timeStamp() @@ -200,10 +209,7 @@ use debug, only: & instance = phase_plasticityInstance(phase) myConstituents = count(material_phase==phase) #ifdef HDF - outID(instance)=HDF5_addGroup(str1,tempResults) -#endif -#ifdef NEWSTATE - allocate(plasticState(phase)%s(1,myConstituents)) ! initialize state (like stateInit) + outID(instance)=HDF5_addGroup(str1,tempResults) #endif endif cycle ! skip to next line @@ -225,18 +231,17 @@ use debug, only: & case ('flowstress') constitutive_j2_outputID(constitutive_j2_Noutput(instance),instance) = flowstress_ID #ifdef HDF - call HDF5_addScalarDataset(outID(instance),a,'flowstress','MPa') - allocate(constitutive_j2_Output2(instance)%flowstress(a)) + call HDF5_addScalarDataset(outID(instance),myConstituents,'flowstress','MPa') + allocate(constitutive_j2_Output2(instance)%flowstress(myConstituents)) constitutive_j2_Output2(instance)%flowstressActive = .true. #endif case ('strainrate') constitutive_j2_outputID(constitutive_j2_Noutput(instance),instance) = strainrate_ID #ifdef HDF - call HDF5_addScalarDataset(outID(instance),a,'strainrate','1/s') - allocate(constitutive_j2_Output2(instance)%strainrate(a)) + call HDF5_addScalarDataset(outID(instance),myConstituents,'strainrate','1/s') + allocate(constitutive_j2_Output2(instance)%strainrate(myConstituents)) constitutive_j2_Output2(instance)%strainrateActive = .true. #endif - case default call IO_error(105_pInt,ext_msg=IO_stringValue(line,positions,2_pInt)//' ('//PLASTICITY_J2_label//')') end select @@ -287,7 +292,9 @@ use debug, only: & enddo parsingFile initializeInstances: do phase = 1_pInt, size(phase_plasticity) - if (phase_plasticity(phase) == PLASTICITY_j2_ID .and. count(material_phase==phase)/=0) then + NofMyPhase=count(material_phase==phase) + if (phase_plasticity(phase) == PLASTICITY_j2_ID .and. NofMyPhase/=0) then + instance = phase_plasticityInstance(phase) outputsLoop: do o = 1_pInt,constitutive_j2_Noutput(instance) select case(constitutive_j2_outputID(o,instance)) case(flowstress_ID,strainrate_ID) @@ -301,15 +308,37 @@ use debug, only: & constitutive_j2_sizePostResults(instance) + mySize endif enddo outputsLoop +#ifdef NEWSTATE + sizeDotState = constitutive_j2_sizeDotState(instance) + sizeState = constitutive_j2_sizeState(instance) + allocate(plasticState(phase)%state0 (sizeState,NofMyPhase),source=constitutive_j2_tau0(instance)) + allocate(plasticState(phase)%partionedState0(sizeState,NofMyPhase),source=constitutive_j2_tau0(instance)) + allocate(plasticState(phase)%subState0 (sizeState,NofMyPhase),source=0.0_pReal) + allocate(plasticState(phase)%state (sizeState,NofMyPhase),source=constitutive_j2_tau0(instance)) + allocate(plasticState(phase)%state_backup (sizeState,NofMyPhase),source=0.0_pReal) + allocate(plasticState(phase)%aTolState (sizeState,NofMyPhase),constitutive_j2_aTolResistance(instance)) + allocate(plasticState(phase)%dotState (sizeDotState,NofMyPhase),source=0.0_pReal) + allocate(plasticState(phase)%deltaState (sizeDotState,NofMyPhase),source=0.0_pReal) + allocate(plasticState(phase)%dotState_backup (sizeDotState,NofMyPhase),source=0.0_pReal) + if (any(numerics_integrator == 1_pInt)) then + allocate(plasticState(phase)%previousDotState (sizeDotState,NofMyPhase),source=0.0_pReal) + allocate(plasticState(phase)%previousDotState2 (sizeDotState,NofMyPhase),source=0.0_pReal) + endif + if (any(numerics_integrator == 4_pInt)) & + allocate(plasticState(phase)%RK4dotState (sizeDotState,NofMyPhase),source=0.0_pReal) + if (any(numerics_integrator == 5_pInt)) & + allocate(plasticState(phase)%RKCK45dotState (sizeDotState,NofMyPhase),source=0.0_pReal) +#endif endif enddo initializeInstances end subroutine constitutive_j2_init -! probably not needed for new state + !-------------------------------------------------------------------------------------------------- !> @brief sets the initial microstructural state for a given instance of this plasticity !> @details initial microstructural state is set to the value specified by tau0 +! not needed for new state !-------------------------------------------------------------------------------------------------- pure function constitutive_j2_stateInit(instance) @@ -325,6 +354,7 @@ end function constitutive_j2_stateInit !-------------------------------------------------------------------------------------------------- !> @brief sets the relevant state values for a given instance of this plasticity +! not needed for new state !-------------------------------------------------------------------------------------------------- pure function constitutive_j2_aTolState(instance) diff --git a/code/constitutive_phenopowerlaw.f90 b/code/constitutive_phenopowerlaw.f90 index 20c9795c1..91897e758 100644 --- a/code/constitutive_phenopowerlaw.f90 +++ b/code/constitutive_phenopowerlaw.f90 @@ -480,7 +480,7 @@ subroutine constitutive_phenopowerlaw_init(fileUnit) constitutive_phenopowerlaw_totalNtwin(instance) ! s_slip, s_twin, sum(gamma), sum(f), accshear_slip, accshear_twin constitutive_phenopowerlaw_sizeState(instance) = constitutive_phenopowerlaw_sizeDotState(instance) - do f = 1_pInt,lattice_maxNslipFamily ! >>> interaction slip -- X + do f = 1_pInt,lattice_maxNslipFamily ! >>> interaction slip -- X index_myFamily = sum(constitutive_phenopowerlaw_Nslip(1:f-1_pInt,instance)) do j = 1_pInt,constitutive_phenopowerlaw_Nslip(f,instance) ! loop over (active) systems in my family (slip) do o = 1_pInt,lattice_maxNslipFamily @@ -505,7 +505,7 @@ subroutine constitutive_phenopowerlaw_init(fileUnit) enddo; enddo - do f = 1_pInt,lattice_maxNtwinFamily ! >>> interaction twin -- X + do f = 1_pInt,lattice_maxNtwinFamily ! >>> interaction twin -- X index_myFamily = sum(constitutive_phenopowerlaw_Ntwin(1:f-1_pInt,instance)) do j = 1_pInt,constitutive_phenopowerlaw_Ntwin(f,instance) ! loop over (active) systems in my family (twin) diff --git a/code/material.f90 b/code/material.f90 index 10fd1ec17..e8e2bd9e3 100644 --- a/code/material.f90 +++ b/code/material.f90 @@ -90,13 +90,14 @@ module material microstructure_crystallite !< crystallite setting ID of each microstructure integer(pInt), dimension(:,:,:), allocatable, public :: & - material_phase !< phase (index) of each grain,IP,element + material_phase !< phase (index) of each grain,IP,element #ifdef NEWSTATE - type(tState), allocatable, dimension(:) :: & + type(tState), allocatable, dimension(:), public :: & plasticState, & elasticState #endif + integer(pInt), dimension(:,:,:), allocatable, public, protected :: & material_texture !< texture (index) of each grain,IP,element @@ -219,10 +220,9 @@ subroutine material_init close(FILEUNIT) #ifdef NEWSTATE allocate(plasticState(material_Nphase)) - allocate(plasticState(material_Nphase)) + allocate(elasticState(material_Nphase)) #endif - - + do m = 1_pInt,material_Nmicrostructure if(microstructure_crystallite(m) < 1_pInt .or. & microstructure_crystallite(m) > material_Ncrystallite) & @@ -845,7 +845,7 @@ subroutine material_populateGrains phaseID,textureID,dGrains,myNgrains,myNorientations,myNconstituents, & grain,constituentGrain,ipGrain,symExtension, ip real(pReal) :: extreme,rnd - integer(pInt), dimension (:,:), allocatable :: Nelems ! counts number of elements in homog, micro array + integer(pInt), dimension (:,:), allocatable :: Nelems ! counts number of elements in homog, micro array type(p_intvec), dimension (:,:), allocatable :: elemsOfHomogMicro ! lists element number in homog, micro array myDebug = debug_level(debug_material) @@ -1157,6 +1157,7 @@ subroutine material_populateGrains end subroutine material_populateGrains +#ifdef HDF integer(pInt) pure function material_NconstituentsPhase(matID) implicit none @@ -1164,5 +1165,6 @@ integer(pInt) pure function material_NconstituentsPhase(matID) material_NconstituentsPhase = count(microstructure_phase == matID) end function +#endif end module material diff --git a/code/prec.f90 b/code/prec.f90 index 6e2f752ca..e68dea7de 100644 --- a/code/prec.f90 +++ b/code/prec.f90 @@ -5,6 +5,7 @@ !> @author Philip Eisenlohr, Max-Planck-Institut für Eisenforschung GmbH !> @author Christoph Kords, Max-Planck-Institut für Eisenforschung GmbH !> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH +!> @author Luv Sharma, Max-Planck-Institut für Eisenforschung GmbH !> @brief setting precision for real and int type depending on makros "FLOAT" and "INT" !> @details setting precision for real and int type and for DAMASK_NaN. Definition is made !! depending on makros "FLOAT" and "INT" defined during compilation @@ -59,10 +60,23 @@ module prec #ifdef NEWSTATE !http://stackoverflow.com/questions/3948210/can-i-have-a-pointer-to-an-item-in-an-allocatable-array type, public :: tState - real(pReal), pointer, dimension(:,:) :: s ! material points, state size + real(pReal), pointer, dimension(:,:) :: state, & ! material points, state size + dotState, & + state0, & + partionedState0, & + subState0, & + state_backup, & + deltaState, & + previousDotState, & + previousDotState2, & + dotState_backup, & + RK4dotState, & + aTolState, & + RKCK45dotState end type #endif + public :: & prec_init