From bd67d2bb6a4c9dd62ad24e6f95f93e1a9384faf7 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 23 Jan 2020 12:48:18 +0100 Subject: [PATCH] new mappings have clear name, not (1,2) --- src/damage_local.f90 | 4 ++-- src/damage_nonlocal.f90 | 2 +- src/homogenization.f90 | 36 ++++++++++++++++----------------- src/homogenization_mech_RGC.f90 | 4 ++-- src/material.f90 | 21 ++----------------- src/thermal_adiabatic.f90 | 4 ++-- src/thermal_conduction.f90 | 4 ++-- 7 files changed, 29 insertions(+), 46 deletions(-) diff --git a/src/damage_local.f90 b/src/damage_local.f90 index 6cb45a391..6296631e5 100644 --- a/src/damage_local.f90 +++ b/src/damage_local.f90 @@ -74,7 +74,7 @@ subroutine damage_local_init allocate(damageState(h)%state (1,NofMyHomog), source=damage_initialPhi(h)) nullify(damageMapping(h)%p) - damageMapping(h)%p => mappingHomogenization(1,:,:) + damageMapping(h)%p => material_homogenizationMemberAt deallocate(damage(h)%p) damage(h)%p => damageState(h)%state(1,:) @@ -103,7 +103,7 @@ function damage_local_updateState(subdt, ip, el) phi, phiDot, dPhiDot_dPhi homog = material_homogenizationAt(el) - offset = mappingHomogenization(1,ip,el) + offset = material_homogenizationMemberAt(ip,el) phi = damageState(homog)%subState0(1,offset) call damage_local_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, el) phi = max(residualStiffness,min(1.0_pReal,phi + subdt*phiDot)) diff --git a/src/damage_nonlocal.f90 b/src/damage_nonlocal.f90 index 17bdecaca..4587857aa 100644 --- a/src/damage_nonlocal.f90 +++ b/src/damage_nonlocal.f90 @@ -79,7 +79,7 @@ subroutine damage_nonlocal_init allocate(damageState(h)%state (1,NofMyHomog), source=damage_initialPhi(h)) nullify(damageMapping(h)%p) - damageMapping(h)%p => mappingHomogenization(1,:,:) + damageMapping(h)%p => material_homogenizationMemberAt deallocate(damage(h)%p) damage(h)%p => damageState(h)%state(1,:) diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 18148e1e0..ba98f87f1 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -242,16 +242,16 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) materialpoint_requested(i,e) = .true. ! everybody requires calculation if (homogState(material_homogenizationAt(e))%sizeState > 0) & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state + homogState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) = & + homogState(material_homogenizationAt(e))%State0( :,material_homogenizationMemberAt(i,e)) ! ...internal homogenization state if (thermalState(material_homogenizationAt(e))%sizeState > 0) & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal thermal state + thermalState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) = & + thermalState(material_homogenizationAt(e))%State0( :,material_homogenizationMemberAt(i,e)) ! ...internal thermal state if (damageState(material_homogenizationAt(e))%sizeState > 0) & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal damage state + damageState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) = & + damageState(material_homogenizationAt(e))%State0( :,material_homogenizationMemberAt(i,e)) ! ...internal damage state enddo enddo @@ -313,14 +313,14 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) enddo if(homogState(material_homogenizationAt(e))%sizeState > 0) & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + homogState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) = & + homogState(material_homogenizationAt(e))%State (:,material_homogenizationMemberAt(i,e)) if(thermalState(material_homogenizationAt(e))%sizeState > 0) & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + thermalState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) = & + thermalState(material_homogenizationAt(e))%State (:,material_homogenizationMemberAt(i,e)) if(damageState(material_homogenizationAt(e))%sizeState > 0) & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%State (:,mappingHomogenization(1,i,e)) + damageState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) = & + damageState(material_homogenizationAt(e))%State (:,material_homogenizationMemberAt(i,e)) materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) @@ -375,14 +375,14 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt) enddo enddo if(homogState(material_homogenizationAt(e))%sizeState > 0) & - homogState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + homogState(material_homogenizationAt(e))%State( :,material_homogenizationMemberAt(i,e)) = & + homogState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) if(thermalState(material_homogenizationAt(e))%sizeState > 0) & - thermalState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + thermalState(material_homogenizationAt(e))%State( :,material_homogenizationMemberAt(i,e)) = & + thermalState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) if(damageState(material_homogenizationAt(e))%sizeState > 0) & - damageState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = & - damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) + damageState(material_homogenizationAt(e))%State( :,material_homogenizationMemberAt(i,e)) = & + damageState(material_homogenizationAt(e))%subState0(:,material_homogenizationMemberAt(i,e)) endif endif converged diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index c493c4190..ec9b94df9 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -108,7 +108,7 @@ module subroutine mech_RGC_init #ifdef DEBUG if (h==material_homogenizationAt(debug_e)) then - prm%of_debug = mappingHomogenization(1,debug_i,debug_e) + prm%of_debug = material_homogenizationMemberAt(debug_i,debug_e) endif #endif @@ -261,7 +261,7 @@ module procedure mech_RGC_updateState endif zeroTimeStep instance = homogenization_typeInstance(material_homogenizationAt(el)) - of = mappingHomogenization(1,ip,el) + of = material_homogenizationMemberAt(ip,el) associate(stt => state(instance), st0 => state0(instance), dst => dependentState(instance), prm => param(instance)) diff --git a/src/material.f90 b/src/material.f90 index a4494ed6e..1d756f7f6 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -113,11 +113,9 @@ module material phase_Nsources, & !< number of source mechanisms active in each phase phase_Nkinematics, & !< number of kinematic mechanisms active in each phase phase_NstiffnessDegradations, & !< number of stiffness degradation mechanisms active in each phase - phase_Noutput, & !< number of '(output)' items per phase phase_elasticityInstance, & !< instance of particular elasticity of each phase phase_plasticityInstance, & !< instance of particular plasticity of each phase homogenization_Ngrains, & !< number of grains in each homogenization - homogenization_Noutput, & !< number of '(output)' items per homogenization homogenization_typeInstance, & !< instance of particular type of each homogenization thermal_typeInstance, & !< instance of particular type of each thermal transport damage_typeInstance !< instance of particular type of each nonlocal damage @@ -129,7 +127,7 @@ module material ! NEW MAPPINGS integer, dimension(:), allocatable, public, protected :: & ! (elem) material_homogenizationAt !< homogenization ID of each element (copy of discretization_homogenizationAt) - integer, dimension(:,:), allocatable, public, protected :: & ! (ip,elem) + integer, dimension(:,:), allocatable, public, target :: & ! (ip,elem) ToDo: ugly target for mapping hack material_homogenizationMemberAt !< position of the element within its homogenization instance integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem) material_phaseAt !< phase ID of each element @@ -177,7 +175,6 @@ module material homogenization_active ! BEGIN DEPRECATED - integer, dimension(:,:,:), allocatable, public, target :: mappingHomogenization !< mapping from material points to offset in heterogenous state/field integer, dimension(:,:), allocatable, private, target :: mappingHomogenizationConst !< mapping from material points to offset in constant state/field ! END DEPRECATED @@ -362,18 +359,8 @@ subroutine material_init !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! BEGIN DEPRECATED - allocate(mappingHomogenization (2,discretization_nIP,discretization_nElem),source=0) allocate(mappingHomogenizationConst( discretization_nIP,discretization_nElem),source=1) - - CounterHomogenization=0 - do e = 1,discretization_nElem - myHomog = discretization_homogenizationAt(e) - do i = 1, discretization_nIP - CounterHomogenization(myHomog) = CounterHomogenization(myHomog) + 1 - mappingHomogenization(1:2,i,e) = [CounterHomogenization(myHomog),huge(1)] - enddo - enddo - + ! hack needed to initialize field values used during constitutive and crystallite initializations do myHomog = 1,size(config_homogenization) thermalMapping (myHomog)%p => mappingHomogenizationConst @@ -401,7 +388,6 @@ subroutine material_parseHomogenization allocate(thermal_typeInstance(size(config_homogenization)), source=0) allocate(damage_typeInstance(size(config_homogenization)), source=0) allocate(homogenization_Ngrains(size(config_homogenization)), source=0) - allocate(homogenization_Noutput(size(config_homogenization)), source=0) allocate(homogenization_active(size(config_homogenization)), source=.false.) !!!!!!!!!!!!!!! allocate(thermal_initialT(size(config_homogenization)), source=300.0_pReal) allocate(damage_initialPhi(size(config_homogenization)), source=1.0_pReal) @@ -411,7 +397,6 @@ subroutine material_parseHomogenization do h=1, size(config_homogenization) - homogenization_Noutput(h) = config_homogenization(h)%countKeys('(output)') tag = config_homogenization(h)%getString('mech') select case (trim(tag)) @@ -548,11 +533,9 @@ subroutine material_parsePhase allocate(phase_Nsources(size(config_phase)), source=0) allocate(phase_Nkinematics(size(config_phase)), source=0) allocate(phase_NstiffnessDegradations(size(config_phase)),source=0) - allocate(phase_Noutput(size(config_phase)), source=0) allocate(phase_localPlasticity(size(config_phase)), source=.false.) do p=1, size(config_phase) - phase_Noutput(p) = config_phase(p)%countKeys('(output)') phase_Nsources(p) = config_phase(p)%countKeys('(source)') phase_Nkinematics(p) = config_phase(p)%countKeys('(kinematics)') phase_NstiffnessDegradations(p) = config_phase(p)%countKeys('(stiffness_degradation)') diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index d96604e59..c6c18a60e 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -77,7 +77,7 @@ subroutine thermal_adiabatic_init allocate(thermalState(h)%state (1,NofMyHomog), source=thermal_initialT(h)) nullify(thermalMapping(h)%p) - thermalMapping(h)%p => mappingHomogenization(1,:,:) + thermalMapping(h)%p => material_homogenizationMemberAt deallocate(temperature(h)%p) temperature(h)%p => thermalState(h)%state(1,:) deallocate(temperatureRate(h)%p) @@ -109,7 +109,7 @@ function thermal_adiabatic_updateState(subdt, ip, el) T, Tdot, dTdot_dT homog = material_homogenizationAt(el) - offset = mappingHomogenization(1,ip,el) + offset = material_homogenizationMemberAt(ip,el) T = thermalState(homog)%subState0(1,offset) call thermal_adiabatic_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index 537b0366b..b182f21ac 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -79,7 +79,7 @@ subroutine thermal_conduction_init allocate(thermalState(h)%state (0,NofMyHomog)) nullify(thermalMapping(h)%p) - thermalMapping(h)%p => mappingHomogenization(1,:,:) + thermalMapping(h)%p => material_homogenizationMemberAt deallocate(temperature (h)%p) allocate (temperature (h)%p(NofMyHomog), source=thermal_initialT(h)) deallocate(temperatureRate(h)%p) @@ -114,7 +114,7 @@ subroutine thermal_conduction_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) constituent homog = material_homogenizationAt(el) - offset = mappingHomogenization(1,ip,el) + offset = material_homogenizationMemberAt(ip,el) instance = thermal_typeInstance(homog) Tdot = 0.0_pReal