From 3dbc07a26bbcc026d84687c3c68760e199ebb158 Mon Sep 17 00:00:00 2001 From: Sharan Roongta Date: Fri, 23 Oct 2020 18:43:52 +0200 Subject: [PATCH 01/36] better error description --- src/IO.f90 | 2 +- src/YAML_types.f90 | 8 ++++---- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/src/IO.f90 b/src/IO.f90 index a95938a7f..260bd94b8 100644 --- a/src/IO.f90 +++ b/src/IO.f90 @@ -491,7 +491,7 @@ subroutine IO_error(error_ID,el,ip,g,instance,ext_msg) case (705) msg = 'Unsupported feature' case (706) - msg = 'Access by incorrect node type' + msg = 'Type mismatch in YAML data node' case (707) msg = 'Abrupt end of file' case (708) diff --git a/src/YAML_types.f90 b/src/YAML_types.f90 index 09a2d0592..f11faf54b 100644 --- a/src/YAML_types.f90 +++ b/src/YAML_types.f90 @@ -311,7 +311,7 @@ function tNode_asScalar(self) result(scalar) class is(tScalar) scalar => self class default - call IO_error(706,ext_msg='tNode_asScalar') + call IO_error(706,ext_msg='Expected "scalar"') end select end function tNode_asScalar @@ -329,7 +329,7 @@ function tNode_asList(self) result(list) class is(tList) list => self class default - call IO_error(706,ext_msg='tNode_asList') + call IO_error(706,ext_msg='Expected "list"') end select end function tNode_asList @@ -347,7 +347,7 @@ function tNode_asDict(self) result(dict) class is(tDict) dict => self class default - call IO_error(706,ext_msg='tNode_asDict') + call IO_error(706,ext_msg='Expected "dict"') end select end function tNode_asDict @@ -641,7 +641,7 @@ function tNode_contains(self,k) result(exists) endif enddo else - call IO_error(706,ext_msg='tNode_contains') + call IO_error(706,ext_msg='Expected "list" or "dict"') endif end function tNode_contains From 68017e49b2e143c2a8c850600412569ebd42ee56 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 12:45:21 +0200 Subject: [PATCH 02/36] consistent name --- src/constitutive.f90 | 4 ++-- src/constitutive_damage.f90 | 2 +- src/constitutive_plastic_nonlocal.f90 | 4 ++-- src/constitutive_thermal.f90 | 2 +- src/crystallite.f90 | 22 ++++++++++---------- src/damage_local.f90 | 4 ++-- src/damage_nonlocal.f90 | 12 +++++------ src/homogenization.f90 | 30 +++++++++++++-------------- src/homogenization_mech_RGC.f90 | 2 +- src/homogenization_mech_isostrain.f90 | 2 +- src/homogenization_mech_none.f90 | 2 +- src/material.f90 | 12 +++++------ src/thermal_adiabatic.f90 | 12 +++++------ src/thermal_conduction.f90 | 16 +++++++------- 14 files changed, 63 insertions(+), 63 deletions(-) diff --git a/src/constitutive.f90 b/src/constitutive.f90 index 64d3f2d31..9bbb9d285 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -127,7 +127,7 @@ module constitutive instance,of,ip,el) real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress - real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem), intent(in) :: & F, & !< deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -753,7 +753,7 @@ function constitutive_collectDotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el of real(pReal), intent(in) :: & subdt !< timestep - real(pReal), intent(in), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: & + real(pReal), intent(in), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem) :: & FArray, & !< elastic deformation gradient FpArray !< plastic deformation gradient real(pReal), intent(in), dimension(3,3) :: & diff --git a/src/constitutive_damage.f90 b/src/constitutive_damage.f90 index 5572192f7..39b8534fe 100644 --- a/src/constitutive_damage.f90 +++ b/src/constitutive_damage.f90 @@ -184,7 +184,7 @@ module subroutine constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi phiDot = 0.0_pReal dPhiDot_dPhi = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) phase = material_phaseAt(grain,el) constituent = material_phasememberAt(grain,ip,el) do source = 1, phase_Nsources(phase) diff --git a/src/constitutive_plastic_nonlocal.f90 b/src/constitutive_plastic_nonlocal.f90 index 8238f17c9..86b7df859 100644 --- a/src/constitutive_plastic_nonlocal.f90 +++ b/src/constitutive_plastic_nonlocal.f90 @@ -976,7 +976,7 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress - real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem), intent(in) :: & F, & !< elastic deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -1176,7 +1176,7 @@ end subroutine plastic_nonlocal_dotState !--------------------------------------------------------------------------------------------------- function rhoDotFlux(F,Fp,timestep, instance,of,ip,el) - real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem), intent(in) :: & F, & !< elastic deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & diff --git a/src/constitutive_thermal.f90 b/src/constitutive_thermal.f90 index 3aefb99a7..944e5d72d 100644 --- a/src/constitutive_thermal.f90 +++ b/src/constitutive_thermal.f90 @@ -95,7 +95,7 @@ module subroutine constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, homog = material_homogenizationAt(el) instance = thermal_typeInstance(homog) - do grain = 1, homogenization_Ngrains(homog) + do grain = 1, homogenization_Nconstituent(homog) phase = material_phaseAt(grain,el) constituent = material_phasememberAt(grain,ip,el) do source = 1, phase_Nsources(phase) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index dc6b9f7da..1d850abaf 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -162,7 +162,7 @@ subroutine crystallite_init debugCrystallite%ip = config_debug%get_asInt('integrationpoint', defaultVal=1) debugCrystallite%grain = config_debug%get_asInt('grain', defaultVal=1) - cMax = homogenization_maxNgrains + cMax = homogenization_maxNconstituent iMax = discretization_nIP eMax = discretization_nElem @@ -253,7 +253,7 @@ subroutine crystallite_init ! initialize !$OMP PARALLEL DO PRIVATE(i,c) do e = FEsolving_execElem(1),FEsolving_execElem(2) - do i = FEsolving_execIP(1), FEsolving_execIP(2); do c = 1, homogenization_Ngrains(material_homogenizationAt(e)) + do i = FEsolving_execIP(1), FEsolving_execIP(2); do c = 1, homogenization_Nconstituent(material_homogenizationAt(e)) crystallite_Fp0(1:3,1:3,c,i,e) = material_orientation0(c,i,e)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005) crystallite_Fp0(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) & / math_det33(crystallite_Fp0(1:3,1:3,c,i,e))**(1.0_pReal/3.0_pReal) @@ -279,7 +279,7 @@ subroutine crystallite_init !$OMP PARALLEL DO do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) call constitutive_dependentState(crystallite_partitionedF0(1:3,1:3,c,i,e), & crystallite_partitionedFp0(1:3,1:3,c,i,e), & c,i,e) ! update dependent state variables to be consistent with basic states @@ -317,7 +317,7 @@ function crystallite_stress() e, & !< counter in element loop startIP, endIP, & s - logical, dimension(homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: todo !ToDo: need to set some values to false for different Ngrains + logical, dimension(homogenization_maxNconstituent,discretization_nIP,discretization_nElem) :: todo !ToDo: need to set some values to false for different Ngrains real(pReal), dimension(:,:,:,:,:), allocatable :: & subLp0,& !< plastic velocity grad at start of crystallite inc subLi0 !< intermediate velocity grad at start of crystallite inc @@ -335,7 +335,7 @@ function crystallite_stress() crystallite_subStep = 0.0_pReal !$OMP PARALLEL DO elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) - do i = FEsolving_execIP(1),FEsolving_execIP(2); do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do i = FEsolving_execIP(1),FEsolving_execIP(2); do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then plasticState (material_phaseAt(c,e))%subState0( :,material_phaseMemberAt(c,i,e)) = & plasticState (material_phaseAt(c,e))%partitionedState0(:,material_phaseMemberAt(c,i,e)) @@ -376,7 +376,7 @@ function crystallite_stress() !$OMP PARALLEL DO PRIVATE(formerSubStep) elementLooping3: do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) !-------------------------------------------------------------------------------------------------- ! wind forward if (crystallite_converged(c,i,e)) then @@ -472,7 +472,7 @@ subroutine crystallite_initializeRestorationPoints(i,e) c, & !< constituent number s - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) crystallite_partitionedFp0(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) crystallite_partitionedLp0(1:3,1:3,c,i,e) = crystallite_Lp0(1:3,1:3,c,i,e) crystallite_partitionedFi0(1:3,1:3,c,i,e) = crystallite_Fi0(1:3,1:3,c,i,e) @@ -503,7 +503,7 @@ subroutine crystallite_windForward(i,e) c, & !< constituent number s - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) crystallite_partitionedF0 (1:3,1:3,c,i,e) = crystallite_partitionedF(1:3,1:3,c,i,e) crystallite_partitionedFp0(1:3,1:3,c,i,e) = crystallite_Fp (1:3,1:3,c,i,e) crystallite_partitionedLp0(1:3,1:3,c,i,e) = crystallite_Lp (1:3,1:3,c,i,e) @@ -536,7 +536,7 @@ subroutine crystallite_restore(i,e,includeL) c, & !< constituent number s - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) if (includeL) then crystallite_Lp(1:3,1:3,c,i,e) = crystallite_partitionedLp0(1:3,1:3,c,i,e) crystallite_Li(1:3,1:3,c,i,e) = crystallite_partitionedLi0(1:3,1:3,c,i,e) @@ -697,7 +697,7 @@ subroutine crystallite_orientations !$OMP PARALLEL DO do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) call crystallite_orientation(c,i,e)%fromMatrix(transpose(math_rotationalPart(crystallite_Fe(1:3,1:3,c,i,e)))) enddo; enddo; enddo !$OMP END PARALLEL DO @@ -848,7 +848,7 @@ subroutine crystallite_results type(rotation), allocatable, dimension(:) :: select_rotations integer :: e,i,c,j - allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNgrains*discretization_nIP)) + allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNconstituent*discretization_nIP)) j=0 do e = 1, size(material_phaseAt,2) diff --git a/src/damage_local.f90 b/src/damage_local.f90 index af2532184..d36139a99 100644 --- a/src/damage_local.f90 +++ b/src/damage_local.f90 @@ -143,8 +143,8 @@ subroutine damage_local_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, el call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) - phiDot = phiDot/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) - dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + phiDot = phiDot/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end subroutine damage_local_getSourceAndItsTangent diff --git a/src/damage_nonlocal.f90 b/src/damage_nonlocal.f90 index f31b10242..f68d778aa 100644 --- a/src/damage_nonlocal.f90 +++ b/src/damage_nonlocal.f90 @@ -110,8 +110,8 @@ subroutine damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, dPhiDot_dPhi = 0.0_pReal call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) - phiDot = phiDot/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) - dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + phiDot = phiDot/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end subroutine damage_nonlocal_getSourceAndItsTangent @@ -132,13 +132,13 @@ function damage_nonlocal_getDiffusion(ip,el) homog = material_homogenizationAt(el) damage_nonlocal_getDiffusion = 0.0_pReal - do grain = 1, homogenization_Ngrains(homog) + do grain = 1, homogenization_Nconstituent(homog) damage_nonlocal_getDiffusion = damage_nonlocal_getDiffusion + & crystallite_push33ToRef(grain,ip,el,lattice_D(1:3,1:3,material_phaseAt(grain,el))) enddo damage_nonlocal_getDiffusion = & - num%charLength**2*damage_nonlocal_getDiffusion/real(homogenization_Ngrains(homog),pReal) + num%charLength**2*damage_nonlocal_getDiffusion/real(homogenization_Nconstituent(homog),pReal) end function damage_nonlocal_getDiffusion @@ -156,12 +156,12 @@ real(pReal) function damage_nonlocal_getMobility(ip,el) damage_nonlocal_getMobility = 0.0_pReal - do ipc = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do ipc = 1, homogenization_Nconstituent(material_homogenizationAt(el)) damage_nonlocal_getMobility = damage_nonlocal_getMobility + lattice_M(material_phaseAt(ipc,el)) enddo damage_nonlocal_getMobility = damage_nonlocal_getMobility/& - real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end function damage_nonlocal_getMobility diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 5aa10fcf3..75f93dc08 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -158,7 +158,7 @@ subroutine homogenization_init debugHomog%grain = config_debug%get_asInt('grain',defaultVal = 1) if (debugHomog%grain < 1 & - .or. debugHomog%grain > homogenization_Ngrains(material_homogenizationAt(debugHomog%element))) & + .or. debugHomog%grain > homogenization_Nconstituent(material_homogenizationAt(debugHomog%element))) & call IO_error(602,ext_msg='constituent', el=debugHomog%element, g=debugHomog%grain) @@ -257,7 +257,7 @@ subroutine materialpoint_stressAndItsTangent(dt) !$OMP PARALLEL DO elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) - myNgrains = homogenization_Ngrains(material_homogenizationAt(e)) + myNgrains = homogenization_Nconstituent(material_homogenizationAt(e)) IpLooping1: do i = FEsolving_execIP(1),FEsolving_execIP(2) if (converged(i,e)) then @@ -327,7 +327,7 @@ subroutine materialpoint_stressAndItsTangent(dt) ! deformation partitioning !$OMP PARALLEL DO PRIVATE(myNgrains) elementLooping2: do e = FEsolving_execElem(1),FEsolving_execElem(2) - myNgrains = homogenization_Ngrains(material_homogenizationAt(e)) + myNgrains = homogenization_Nconstituent(material_homogenizationAt(e)) IpLooping2: do i = FEsolving_execIP(1),FEsolving_execIP(2) if(requested(i,e) .and. .not. doneAndHappy(1,i,e)) then ! requested but not yet done call partitionDeformation(materialpoint_F0(1:3,1:3,i,e) & @@ -408,12 +408,12 @@ subroutine partitionDeformation(subF,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization call mech_isostrain_partitionDeformation(& - crystallite_partitionedF(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & subF) case (HOMOGENIZATION_RGC_ID) chosenHomogenization call mech_RGC_partitionDeformation(& - crystallite_partitionedF(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & subF,& ip, & el) @@ -437,19 +437,19 @@ function updateState(subdt,subF,ip,el) el !< element number integer :: c logical, dimension(2) :: updateState - real(pReal) :: dPdFs(3,3,3,3,homogenization_Ngrains(material_homogenizationAt(el))) + real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituent(material_homogenizationAt(el))) updateState = .true. chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - do c=1,homogenization_Ngrains(material_homogenizationAt(el)) + do c=1,homogenization_Nconstituent(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo updateState = & updateState .and. & - mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & - crystallite_partitionedF(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & - crystallite_partitionedF0(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el),& + mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF0(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el),& subF,& subdt, & dPdFs, & @@ -487,7 +487,7 @@ subroutine averageStressAndItsTangent(ip,el) ip, & !< integration point el !< element number integer :: c - real(pReal) :: dPdFs(3,3,3,3,homogenization_Ngrains(material_homogenizationAt(el))) + real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituent(material_homogenizationAt(el))) chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) @@ -496,24 +496,24 @@ subroutine averageStressAndItsTangent(ip,el) materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_stressTangent(1,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization - do c = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do c = 1, homogenization_Nconstituent(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_isostrain_averageStressAndItsTangent(& materialpoint_P(1:3,1:3,ip,el), & materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& - crystallite_P(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - do c = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do c = 1, homogenization_Nconstituent(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_RGC_averageStressAndItsTangent(& materialpoint_P(1:3,1:3,ip,el), & materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& - crystallite_P(1:3,1:3,1:homogenization_Ngrains(material_homogenizationAt(el)),ip,el), & + crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) end select chosenHomogenization diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index ea5bffd28..6bbbbbaa9 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -164,7 +164,7 @@ module subroutine mech_RGC_init(num_homogMech) #endif prm%N_constituents = homogMech%get_asInts('cluster_size',requiredSize=3) - if (homogenization_Ngrains(h) /= product(prm%N_constituents)) & + if (homogenization_Nconstituent(h) /= product(prm%N_constituents)) & call IO_error(211,ext_msg='N_constituents (mech_RGC)') prm%xi_alpha = homogMech%get_asFloat('xi_alpha') diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index f064578c8..5a03e1204 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -51,7 +51,7 @@ module subroutine mech_isostrain_init homogMech => homog%get('mech') associate(prm => param(homogenization_typeInstance(h))) - prm%N_constituents = homogenization_Ngrains(h) + prm%N_constituents = homogenization_Nconstituent(h) select case(homogMech%get_asString('mapping',defaultVal = 'sum')) case ('sum') prm%mapping = parallel_ID diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index a58147c45..d64154dcc 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -26,7 +26,7 @@ module subroutine mech_none_init do h = 1, size(homogenization_type) if(homogenization_type(h) /= HOMOGENIZATION_NONE_ID) cycle - if(homogenization_Ngrains(h) /= 1) & + if(homogenization_Nconstituent(h) /= 1) & call IO_error(211,ext_msg='N_constituents (mech_none)') NofMyHomog = count(material_homogenizationAt == h) diff --git a/src/material.f90 b/src/material.f90 index 11dfeb42e..ae6630f00 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -67,10 +67,10 @@ module material material_Nhomogenization !< number of homogenizations integer, public, protected :: & - homogenization_maxNgrains !< max number of grains in any USED homogenization + homogenization_maxNconstituent !< max number of grains in any USED homogenization integer, dimension(:), allocatable, public, protected :: & - homogenization_Ngrains, & !< number of grains in each homogenization + homogenization_Nconstituent, & !< number of grains in each 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 @@ -187,7 +187,7 @@ subroutine material_init(restart) print*, 'Homogenization parsed' - if(homogenization_maxNgrains > size(material_phaseAt,1)) call IO_error(148) + if(homogenization_maxNconstituent > size(material_phaseAt,1)) call IO_error(148) allocate(homogState (material_Nhomogenization)) allocate(thermalState (material_Nhomogenization)) @@ -249,14 +249,14 @@ subroutine material_parseHomogenization allocate(homogenization_typeInstance(material_Nhomogenization), source=0) allocate(thermal_typeInstance(material_Nhomogenization), source=0) allocate(damage_typeInstance(material_Nhomogenization), source=0) - allocate(homogenization_Ngrains(material_Nhomogenization), source=0) + allocate(homogenization_Nconstituent(material_Nhomogenization), source=0) allocate(thermal_initialT(material_Nhomogenization), source=300.0_pReal) allocate(damage_initialPhi(material_Nhomogenization), source=1.0_pReal) do h=1, material_Nhomogenization homog => material_homogenization%get(h) homogMech => homog%get('mech') - homogenization_Ngrains(h) = homog%get_asInt('N_constituents') + homogenization_Nconstituent(h) = homog%get_asInt('N_constituents') select case (homogMech%get_asString('type')) case('none') homogenization_type(h) = HOMOGENIZATION_NONE_ID @@ -308,7 +308,7 @@ subroutine material_parseHomogenization damage_typeInstance(h) = count(damage_type (1:h) == damage_type (h)) enddo - homogenization_maxNgrains = maxval(homogenization_Ngrains) + homogenization_maxNconstituent = maxval(homogenization_Nconstituent) end subroutine material_parseHomogenization diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index 63deb3cd5..00b73862a 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -145,8 +145,8 @@ subroutine thermal_adiabatic_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) homog = material_homogenizationAt(el) call constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, crystallite_S, crystallite_Lp, ip, el) - Tdot = Tdot/real(homogenization_Ngrains(homog),pReal) - dTdot_dT = dTdot_dT/real(homogenization_Ngrains(homog),pReal) + Tdot = Tdot/real(homogenization_Nconstituent(homog),pReal) + dTdot_dT = dTdot_dT/real(homogenization_Nconstituent(homog),pReal) end subroutine thermal_adiabatic_getSourceAndItsTangent @@ -167,13 +167,13 @@ function thermal_adiabatic_getSpecificHeat(ip,el) thermal_adiabatic_getSpecificHeat = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) thermal_adiabatic_getSpecificHeat = thermal_adiabatic_getSpecificHeat & + lattice_c_p(material_phaseAt(grain,el)) enddo thermal_adiabatic_getSpecificHeat = thermal_adiabatic_getSpecificHeat & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end function thermal_adiabatic_getSpecificHeat @@ -193,13 +193,13 @@ function thermal_adiabatic_getMassDensity(ip,el) thermal_adiabatic_getMassDensity = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) thermal_adiabatic_getMassDensity = thermal_adiabatic_getMassDensity & + lattice_rho(material_phaseAt(grain,el)) enddo thermal_adiabatic_getMassDensity = thermal_adiabatic_getMassDensity & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end function thermal_adiabatic_getMassDensity diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index 69ce8025e..f1f444c28 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -104,8 +104,8 @@ subroutine thermal_conduction_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) homog = material_homogenizationAt(el) call constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, crystallite_S,crystallite_Lp ,ip, el) - Tdot = Tdot/real(homogenization_Ngrains(homog),pReal) - dTdot_dT = dTdot_dT/real(homogenization_Ngrains(homog),pReal) + Tdot = Tdot/real(homogenization_Nconstituent(homog),pReal) + dTdot_dT = dTdot_dT/real(homogenization_Nconstituent(homog),pReal) end subroutine thermal_conduction_getSourceAndItsTangent @@ -125,13 +125,13 @@ function thermal_conduction_getConductivity(ip,el) thermal_conduction_getConductivity = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) thermal_conduction_getConductivity = thermal_conduction_getConductivity + & crystallite_push33ToRef(grain,ip,el,lattice_K(:,:,material_phaseAt(grain,el))) enddo thermal_conduction_getConductivity = thermal_conduction_getConductivity & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end function thermal_conduction_getConductivity @@ -151,13 +151,13 @@ function thermal_conduction_getSpecificHeat(ip,el) thermal_conduction_getSpecificHeat = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) thermal_conduction_getSpecificHeat = thermal_conduction_getSpecificHeat & + lattice_c_p(material_phaseAt(grain,el)) enddo thermal_conduction_getSpecificHeat = thermal_conduction_getSpecificHeat & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end function thermal_conduction_getSpecificHeat @@ -178,13 +178,13 @@ function thermal_conduction_getMassDensity(ip,el) thermal_conduction_getMassDensity = 0.0_pReal - do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) thermal_conduction_getMassDensity = thermal_conduction_getMassDensity & + lattice_rho(material_phaseAt(grain,el)) enddo thermal_conduction_getMassDensity = thermal_conduction_getMassDensity & - / real(homogenization_Ngrains(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) end function thermal_conduction_getMassDensity From 5f2512e4d5e5c44820646fb394c8fdbe66654b9b Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 14:23:00 +0200 Subject: [PATCH 03/36] prepare for more sensible separation of functionality --- src/material.f90 | 56 ++++++++++++++++++++++++++++++++++-------------- 1 file changed, 40 insertions(+), 16 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index ae6630f00..891d786a3 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -67,10 +67,10 @@ module material material_Nhomogenization !< number of homogenizations integer, public, protected :: & - homogenization_maxNconstituent !< max number of grains in any USED homogenization + homogenization_maxNconstituent !< max number of grains in any USED homogenization integer, dimension(:), allocatable, public, protected :: & - homogenization_Nconstituent, & !< number of grains in each homogenization + homogenization_Nconstituent, & !< number of grains in each 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 @@ -85,7 +85,7 @@ module material material_homogenizationMemberAt !< position of the element within its homogenization instance integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem) material_phaseAt !< phase ID of each element - integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,elem) + integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,IP,elem) material_phaseMemberAt !< position of the element within its phase instance type(tState), allocatable, dimension(:), public :: & @@ -183,6 +183,7 @@ subroutine material_init(restart) call material_parseMaterial print*, 'Material parsed' + call material_parseNconstituent call material_parseHomogenization print*, 'Homogenization parsed' @@ -225,6 +226,34 @@ subroutine material_init(restart) end subroutine material_init + +!-------------------------------------------------------------------------------------------------- +!> @brief Determine Nconstituent in homogenization +!-------------------------------------------------------------------------------------------------- +subroutine material_parseNconstituent + + class(tNode), pointer :: & + material_homogenization, & + homog + + integer :: h + + material_homogenization => config_material%get('homogenization') + material_Nhomogenization = material_homogenization%length + + allocate(homogenization_Nconstituent(material_Nhomogenization)) + + + do h=1, material_Nhomogenization + homog => material_homogenization%get(h) + homogenization_Nconstituent(h) = homog%get_asInt('N_constituents') + enddo + homogenization_maxNconstituent = maxval(homogenization_Nconstituent) + + +end subroutine material_parseNconstituent + + !-------------------------------------------------------------------------------------------------- !> @brief parses the homogenization part from the material configuration ! ToDo: This should be done in homogenization @@ -241,7 +270,6 @@ subroutine material_parseHomogenization integer :: h material_homogenization => config_material%get('homogenization') - material_Nhomogenization = material_homogenization%length allocate(homogenization_type(material_Nhomogenization), source=HOMOGENIZATION_undefined_ID) allocate(thermal_type(material_Nhomogenization), source=THERMAL_isothermal_ID) @@ -249,7 +277,6 @@ subroutine material_parseHomogenization allocate(homogenization_typeInstance(material_Nhomogenization), source=0) allocate(thermal_typeInstance(material_Nhomogenization), source=0) allocate(damage_typeInstance(material_Nhomogenization), source=0) - allocate(homogenization_Nconstituent(material_Nhomogenization), source=0) allocate(thermal_initialT(material_Nhomogenization), source=300.0_pReal) allocate(damage_initialPhi(material_Nhomogenization), source=1.0_pReal) @@ -308,9 +335,6 @@ subroutine material_parseHomogenization damage_typeInstance(h) = count(damage_type (1:h) == damage_type (h)) enddo - homogenization_maxNconstituent = maxval(homogenization_Nconstituent) - - end subroutine material_parseHomogenization @@ -350,14 +374,14 @@ subroutine material_parseMaterial material_Nconstituents(m) = constituents%length enddo maxNconstituents = maxval(material_Nconstituents) - + allocate(material_homogenizationAt(discretization_nElem),source=0) allocate(material_homogenizationMemberAt(discretization_nIP,discretization_nElem),source=0) allocate(material_phaseAt(maxNconstituents,discretization_nElem),source=0) allocate(material_phaseMemberAt(maxNconstituents,discretization_nIP,discretization_nElem),source=0) allocate(material_orientation0(maxNconstituents,discretization_nIP,discretization_nElem)) - + phases => config_material%get('phase') allocate(counterPhase(phases%length),source=0) homogenizations => config_material%get('homogenization') @@ -366,29 +390,29 @@ subroutine material_parseMaterial do e = 1, discretization_nElem material => materials%get(discretization_materialAt(e)) constituents => material%get('constituents') - + material_homogenizationAt(e) = homogenizations%getIndex(material%get_asString('homogenization')) do i = 1, discretization_nIP counterHomogenization(material_homogenizationAt(e)) = counterHomogenization(material_homogenizationAt(e)) + 1 material_homogenizationMemberAt(i,e) = counterHomogenization(material_homogenizationAt(e)) enddo - + frac = 0.0_pReal do c = 1, constituents%length constituent => constituents%get(c) frac = frac + constituent%get_asFloat('fraction') - + material_phaseAt(c,e) = phases%getIndex(constituent%get_asString('phase')) do i = 1, discretization_nIP counterPhase(material_phaseAt(c,e)) = counterPhase(material_phaseAt(c,e)) + 1 material_phaseMemberAt(c,i,e) = counterPhase(material_phaseAt(c,e)) - + call material_orientation0(c,i,e)%fromQuaternion(constituent%get_asFloats('O',requiredSize=4)) enddo - + enddo if (dNeq(frac,1.0_pReal)) call IO_error(153,ext_msg='constituent') - + enddo end subroutine material_parseMaterial From 5e4815cddfdaf0f1974b05d60302863aae1870c3 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 14:48:07 +0200 Subject: [PATCH 04/36] not needed as public variable --- src/material.f90 | 25 +++++++++++++------------ 1 file changed, 13 insertions(+), 12 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index 891d786a3..f16c10362 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -52,7 +52,7 @@ module material HOMOGENIZATION_RGC_ID end enum - character(len=pStringLen), public, protected, allocatable, dimension(:) :: & + character(len=pStringLen), public, protected, allocatable, dimension(:) :: & material_name_phase, & !< name of each phase material_name_homogenization !< name of each homogenization @@ -96,11 +96,6 @@ module material type(Rotation), dimension(:,:,:), allocatable, public, protected :: & material_orientation0 !< initial orientation of each grain,IP,element - integer, dimension(:), allocatable, private :: & - material_Nconstituents !< number of constituents in each material - - - ! BEGIN DEPRECATED integer, dimension(:,:), allocatable, private, target :: mappingHomogenizationConst !< mapping from material points to offset in constant state/field ! END DEPRECATED @@ -363,14 +358,23 @@ subroutine material_parseMaterial c, & maxNconstituents + integer, dimension(:), allocatable :: & + material_Nconstituents !< number of constituents in each material + materials => config_material%get('material') if(any(discretization_materialAt > materials%length)) & call IO_error(155,ext_msg='More materials requested than found in material.yaml') + phases => config_material%get('phase') + allocate(counterPhase(phases%length),source=0) + homogenizations => config_material%get('homogenization') + allocate(counterHomogenization(homogenizations%length),source=0) + allocate(material_Nconstituents(materials%length),source=0) + do m = 1, materials%length material => materials%get(m) - constituents => material%get('constituents') + constituents => material%get('constituents') material_Nconstituents(m) = constituents%length enddo maxNconstituents = maxval(material_Nconstituents) @@ -382,10 +386,7 @@ subroutine material_parseMaterial allocate(material_orientation0(maxNconstituents,discretization_nIP,discretization_nElem)) - phases => config_material%get('phase') - allocate(counterPhase(phases%length),source=0) - homogenizations => config_material%get('homogenization') - allocate(counterHomogenization(homogenizations%length),source=0) + do e = 1, discretization_nElem material => materials%get(discretization_materialAt(e)) @@ -407,7 +408,7 @@ subroutine material_parseMaterial counterPhase(material_phaseAt(c,e)) = counterPhase(material_phaseAt(c,e)) + 1 material_phaseMemberAt(c,i,e) = counterPhase(material_phaseAt(c,e)) - call material_orientation0(c,i,e)%fromQuaternion(constituent%get_asFloats('O',requiredSize=4)) + call material_orientation0(c,i,e)%fromQuaternion(constituent%get_asFloats('O',requiredSize=4)) ! should be done in crystallite enddo enddo From a1ab526456be14804f20159f0e638f110fde9bb6 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 14:56:03 +0200 Subject: [PATCH 05/36] avoid global variables --- src/crystallite.f90 | 6 ++--- src/damage_none.f90 | 2 +- src/grid/DAMASK_grid.f90 | 2 +- src/material.f90 | 50 +++++++++++++++++--------------------- src/thermal_adiabatic.f90 | 2 +- src/thermal_conduction.f90 | 2 +- src/thermal_isothermal.f90 | 2 +- 7 files changed, 30 insertions(+), 36 deletions(-) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index 1d850abaf..b71bb37f0 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -1571,7 +1571,7 @@ subroutine crystallite_restartWrite call HDF5_closeGroup(groupHandle) groupHandle = HDF5_addGroup(fileHandle,'materialpoint') - do i = 1, material_Nhomogenization + do i = 1, size(material_name_homogenization) write(datasetName,'(i0,a)') i,'_omega_homogenization' call HDF5_write(groupHandle,homogState(i)%state,datasetName) enddo @@ -1612,7 +1612,7 @@ subroutine crystallite_restartRead call HDF5_closeGroup(groupHandle) groupHandle = HDF5_openGroup(fileHandle,'materialpoint') - do i = 1, material_Nhomogenization + do i = 1,size(material_name_homogenization) write(datasetName,'(i0,a)') i,'_omega_homogenization' call HDF5_read(groupHandle,homogState(i)%state0,datasetName) enddo @@ -1645,7 +1645,7 @@ subroutine crystallite_forward do j = 1,phase_Nsources(i) sourceState(i)%p(j)%state0 = sourceState(i)%p(j)%state enddo; enddo - do i = 1, material_Nhomogenization + do i = 1,size(material_name_homogenization) homogState (i)%state0 = homogState (i)%state thermalState(i)%state0 = thermalState(i)%state damageState (i)%state0 = damageState (i)%state diff --git a/src/damage_none.f90 b/src/damage_none.f90 index 4e01998bc..52100707e 100644 --- a/src/damage_none.f90 +++ b/src/damage_none.f90 @@ -20,7 +20,7 @@ subroutine damage_none_init print'(/,a)', ' <<<+- damage_none init -+>>>'; flush(6) - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (damage_type(h) /= DAMAGE_NONE_ID) cycle NofMyHomog = count(material_homogenizationAt == h) diff --git a/src/grid/DAMASK_grid.f90 b/src/grid/DAMASK_grid.f90 index 4b8302def..1cb3acf16 100644 --- a/src/grid/DAMASK_grid.f90 +++ b/src/grid/DAMASK_grid.f90 @@ -210,7 +210,7 @@ program DAMASK_grid if(.not. step_discretization%contains('t')) call IO_error(error_ID=837,ext_msg = 't missing') if(.not. step_discretization%contains('N')) call IO_error(error_ID=837,ext_msg = 'N missing') loadCases(l)%time = step_discretization%get_asFloat('t') - loadCases(l)%incs = step_discretization%get_asFloat('N') + loadCases(l)%incs = step_discretization%get_asInt ('N') loadCases(l)%logscale = step_discretization%get_asBool ('log_timestep', defaultVal= .false.) loadCases(l)%outputfrequency = step_discretization%get_asInt ('f_out', defaultVal=1) loadCases(l)%restartfrequency = step_discretization%get_asInt ('f_restart', defaultVal=huge(0)) diff --git a/src/material.f90 b/src/material.f90 index f16c10362..d65adf7c7 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -63,9 +63,6 @@ module material integer(kind(HOMOGENIZATION_undefined_ID)), dimension(:), allocatable, public, protected :: & homogenization_type !< type of each homogenization - integer, public, protected :: & - material_Nhomogenization !< number of homogenizations - integer, public, protected :: & homogenization_maxNconstituent !< max number of grains in any USED homogenization @@ -83,7 +80,7 @@ module material material_homogenizationAt !< homogenization ID of each element 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) + integer, dimension(:,:), allocatable, public, protected :: & ! (constituent,elem) material_phaseAt !< phase ID of each element integer, dimension(:,:,:), allocatable, public, protected :: & ! (constituent,IP,elem) material_phaseMemberAt !< position of the element within its phase instance @@ -185,17 +182,17 @@ subroutine material_init(restart) if(homogenization_maxNconstituent > size(material_phaseAt,1)) call IO_error(148) - allocate(homogState (material_Nhomogenization)) - allocate(thermalState (material_Nhomogenization)) - allocate(damageState (material_Nhomogenization)) + allocate(homogState (size(material_name_homogenization))) + allocate(thermalState (size(material_name_homogenization))) + allocate(damageState (size(material_name_homogenization))) - allocate(thermalMapping (material_Nhomogenization)) - allocate(damageMapping (material_Nhomogenization)) + allocate(thermalMapping (size(material_name_homogenization))) + allocate(damageMapping (size(material_name_homogenization))) - allocate(temperature (material_Nhomogenization)) - allocate(damage (material_Nhomogenization)) + allocate(temperature (size(material_name_homogenization))) + allocate(damage (size(material_name_homogenization))) - allocate(temperatureRate (material_Nhomogenization)) + allocate(temperatureRate (size(material_name_homogenization))) if (.not. restart) then @@ -210,7 +207,7 @@ subroutine material_init(restart) allocate(mappingHomogenizationConst( discretization_nIP,discretization_nElem),source=1) ! hack needed to initialize field values used during constitutive initialization - do myHomog = 1,material_Nhomogenization + do myHomog = 1, size(material_name_homogenization) thermalMapping (myHomog)%p => mappingHomogenizationConst damageMapping (myHomog)%p => mappingHomogenizationConst allocate(temperature (myHomog)%p(1), source=thermal_initialT(myHomog)) @@ -234,12 +231,11 @@ subroutine material_parseNconstituent integer :: h material_homogenization => config_material%get('homogenization') - material_Nhomogenization = material_homogenization%length - allocate(homogenization_Nconstituent(material_Nhomogenization)) + allocate(homogenization_Nconstituent(size(material_name_homogenization))) - do h=1, material_Nhomogenization + do h=1, size(material_name_homogenization) homog => material_homogenization%get(h) homogenization_Nconstituent(h) = homog%get_asInt('N_constituents') enddo @@ -266,19 +262,18 @@ subroutine material_parseHomogenization material_homogenization => config_material%get('homogenization') - allocate(homogenization_type(material_Nhomogenization), source=HOMOGENIZATION_undefined_ID) - allocate(thermal_type(material_Nhomogenization), source=THERMAL_isothermal_ID) - allocate(damage_type (material_Nhomogenization), source=DAMAGE_none_ID) - allocate(homogenization_typeInstance(material_Nhomogenization), source=0) - allocate(thermal_typeInstance(material_Nhomogenization), source=0) - allocate(damage_typeInstance(material_Nhomogenization), source=0) - allocate(thermal_initialT(material_Nhomogenization), source=300.0_pReal) - allocate(damage_initialPhi(material_Nhomogenization), source=1.0_pReal) + allocate(homogenization_type(size(material_name_homogenization)), source=HOMOGENIZATION_undefined_ID) + allocate(thermal_type(size(material_name_homogenization)), source=THERMAL_isothermal_ID) + allocate(damage_type (size(material_name_homogenization)), source=DAMAGE_none_ID) + allocate(homogenization_typeInstance(size(material_name_homogenization)), source=0) + allocate(thermal_typeInstance(size(material_name_homogenization)), source=0) + allocate(damage_typeInstance(size(material_name_homogenization)), source=0) + allocate(thermal_initialT(size(material_name_homogenization)), source=300.0_pReal) + allocate(damage_initialPhi(size(material_name_homogenization)), source=1.0_pReal) - do h=1, material_Nhomogenization + do h=1, size(material_name_homogenization) homog => material_homogenization%get(h) homogMech => homog%get('mech') - homogenization_Nconstituent(h) = homog%get_asInt('N_constituents') select case (homogMech%get_asString('type')) case('none') homogenization_type(h) = HOMOGENIZATION_NONE_ID @@ -324,7 +319,7 @@ subroutine material_parseHomogenization endif enddo - do h=1, material_Nhomogenization + do h=1, size(material_name_homogenization) homogenization_typeInstance(h) = count(homogenization_type(1:h) == homogenization_type(h)) thermal_typeInstance(h) = count(thermal_type (1:h) == thermal_type (h)) damage_typeInstance(h) = count(damage_type (1:h) == damage_type (h)) @@ -387,7 +382,6 @@ subroutine material_parseMaterial allocate(material_orientation0(maxNconstituents,discretization_nIP,discretization_nElem)) - do e = 1, discretization_nElem material => materials%get(discretization_materialAt(e)) constituents => material%get('constituents') diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index 00b73862a..07dbdebff 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -54,7 +54,7 @@ subroutine thermal_adiabatic_init allocate(param(maxNinstance)) material_homogenization => config_material%get('homogenization') - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_adiabatic_ID) cycle homog => material_homogenization%get(h) homogThermal => homog%get('thermal') diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index f1f444c28..616965df0 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -53,7 +53,7 @@ subroutine thermal_conduction_init allocate(param(Ninstance)) material_homogenization => config_material%get('homogenization') - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_conduction_ID) cycle homog => material_homogenization%get(h) homogThermal => homog%get('thermal') diff --git a/src/thermal_isothermal.f90 b/src/thermal_isothermal.f90 index 2dda358ac..703a9aaac 100644 --- a/src/thermal_isothermal.f90 +++ b/src/thermal_isothermal.f90 @@ -20,7 +20,7 @@ subroutine thermal_isothermal_init print'(/,a)', ' <<<+- thermal_isothermal init -+>>>'; flush(6) - do h = 1, material_Nhomogenization + do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_isothermal_ID) cycle NofMyHomog = count(material_homogenizationAt == h) From d765a17f0a125788bc7147a5148b70af1f41e933 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 16:30:48 +0200 Subject: [PATCH 06/36] use function to parse material.yaml --- src/material.f90 | 111 ++++++++++++++++++----------------------------- 1 file changed, 42 insertions(+), 69 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index d65adf7c7..7197506f4 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -149,39 +149,18 @@ contains subroutine material_init(restart) logical, intent(in) :: restart - - integer :: ph, myHomog - class(tNode), pointer :: & - phases, & - material_homogenization - character(len=pStringLen) :: sectionName + integer :: myHomog print'(/,a)', ' <<<+- material init -+>>>'; flush(IO_STDOUT) - phases => config_material%get('phase') - allocate(material_name_phase(phases%length)) - do ph = 1, phases%length - write(sectionName,'(i0,a)') ph,'_' - material_name_phase(ph) = trim(adjustl(sectionName))//phases%getKey(ph) !ToDO: No reason to do. Update damage tests - enddo - - material_homogenization => config_material%get('homogenization') - allocate(material_name_homogenization(material_homogenization%length)) - do myHomog = 1, material_homogenization%length - write(sectionName,'(i0,a)') myHomog,'_' - material_name_homogenization(myHomog) = trim(adjustl(sectionName))//material_homogenization%getKey(myHomog) - enddo call material_parseMaterial print*, 'Material parsed' - call material_parseNconstituent call material_parseHomogenization print*, 'Homogenization parsed' - if(homogenization_maxNconstituent > size(material_phaseAt,1)) call IO_error(148) - allocate(homogState (size(material_name_homogenization))) allocate(thermalState (size(material_name_homogenization))) allocate(damageState (size(material_name_homogenization))) @@ -219,32 +198,6 @@ subroutine material_init(restart) end subroutine material_init -!-------------------------------------------------------------------------------------------------- -!> @brief Determine Nconstituent in homogenization -!-------------------------------------------------------------------------------------------------- -subroutine material_parseNconstituent - - class(tNode), pointer :: & - material_homogenization, & - homog - - integer :: h - - material_homogenization => config_material%get('homogenization') - - allocate(homogenization_Nconstituent(size(material_name_homogenization))) - - - do h=1, size(material_name_homogenization) - homog => material_homogenization%get(h) - homogenization_Nconstituent(h) = homog%get_asInt('N_constituents') - enddo - homogenization_maxNconstituent = maxval(homogenization_Nconstituent) - - -end subroutine material_parseNconstituent - - !-------------------------------------------------------------------------------------------------- !> @brief parses the homogenization part from the material configuration ! ToDo: This should be done in homogenization @@ -338,49 +291,69 @@ subroutine material_parseMaterial constituents, & !> list of constituents constituent, & !> constituent definition phases, & - homogenizations + homogenizations, & + homogenization integer, dimension(:), allocatable :: & counterPhase, & counterHomogenization + character(len=pStringLen) :: sectionName + real(pReal) :: & frac integer :: & - e, & - i, & - m, & - c, & - maxNconstituents + e, i, c, & + p, h, & + m + + materials => config_material%get('material') + phases => config_material%get('phase') + homogenizations => config_material%get('homogenization') - integer, dimension(:), allocatable :: & - material_Nconstituents !< number of constituents in each material - materials => config_material%get('material') if(any(discretization_materialAt > materials%length)) & call IO_error(155,ext_msg='More materials requested than found in material.yaml') - phases => config_material%get('phase') - allocate(counterPhase(phases%length),source=0) - homogenizations => config_material%get('homogenization') - allocate(counterHomogenization(homogenizations%length),source=0) +!-------------------------------------------------------------------------------------------------- +! determine name of phases/homogenizations + allocate(material_name_phase(phases%length)) + do p = 1, phases%length + write(sectionName,'(i0,a)') p,'_' + material_name_phase(p) = trim(adjustl(sectionName))//phases%getKey(p) !ToDo: remove prefix + enddo - allocate(material_Nconstituents(materials%length),source=0) + allocate(material_name_homogenization(homogenizations%length)) + do h = 1, homogenizations%length + write(sectionName,'(i0,a)') h,'_' + material_name_homogenization(h) = trim(adjustl(sectionName))//homogenizations%getKey(h) !ToDo: remove prefix + enddo +!-------------------------------------------------------------------------------------------------- +! sanity check: Matching # of constituents do m = 1, materials%length material => materials%get(m) - constituents => material%get('constituents') - material_Nconstituents(m) = constituents%length + constituents => material%get('constituents') + homogenization => homogenizations%get(material%get_asString('homogenization')) + if(constituents%length /= homogenization%get_asInt('N_constituents')) call IO_error(148) enddo - maxNconstituents = maxval(material_Nconstituents) + + allocate(homogenization_Nconstituent(size(material_name_homogenization))) + do h=1, size(material_name_homogenization) + homogenization => homogenizations%get(h) + homogenization_Nconstituent(h) = homogenization%get_asInt('N_constituents') + enddo + homogenization_maxNconstituent = maxval(homogenization_Nconstituent) + + allocate(counterPhase(phases%length),source=0) + allocate(counterHomogenization(homogenizations%length),source=0) allocate(material_homogenizationAt(discretization_nElem),source=0) allocate(material_homogenizationMemberAt(discretization_nIP,discretization_nElem),source=0) - allocate(material_phaseAt(maxNconstituents,discretization_nElem),source=0) - allocate(material_phaseMemberAt(maxNconstituents,discretization_nIP,discretization_nElem),source=0) - - allocate(material_orientation0(maxNconstituents,discretization_nIP,discretization_nElem)) + allocate(material_phaseAt(homogenization_maxNconstituent,discretization_nElem),source=0) + allocate(material_phaseMemberAt(homogenization_maxNconstituent,discretization_nIP,discretization_nElem),source=0) + allocate(material_orientation0(homogenization_maxNconstituent,discretization_nIP,discretization_nElem)) do e = 1, discretization_nElem material => materials%get(discretization_materialAt(e)) From e464f11412bbc5d8bd3d7dc548adfac649b59a08 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 17:26:24 +0200 Subject: [PATCH 07/36] polishing --- src/YAML_types.f90 | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/YAML_types.f90 b/src/YAML_types.f90 index 09a2d0592..cd4be5c7e 100644 --- a/src/YAML_types.f90 +++ b/src/YAML_types.f90 @@ -583,9 +583,9 @@ function tNode_get_byIndex_asStrings(self,i) result(nodeAsStrings) end function tNode_get_byIndex_asStrings -!------------------------------------------------------------------------------------------------------- +!-------------------------------------------------------------------------------------------------- !> @brief Returns the key in a dictionary as a string -!------------------------------------------------------------------------------------------------------- +!-------------------------------------------------------------------------------------------------- function tNode_getKey_byIndex(self,i) result(key) class(tNode), intent(in), target :: self @@ -1169,17 +1169,21 @@ function tList_asStrings(self) type(tScalar), pointer :: scalar len_max = 0 - allocate(character(len=pStringLen) :: tList_asStrings(self%length)) + item => self%first + do i = 1, self%length + scalar => item%node%asScalar() + len_max = max(len_max, len_trim(scalar%asString())) + item => item%next + enddo + + allocate(character(len=len_max) :: tList_asStrings(self%length)) item => self%first do i = 1, self%length scalar => item%node%asScalar() tList_asStrings(i) = scalar%asString() - len_max = max(len_max, len_trim(tList_asStrings(i))) item => item%next enddo - !ToDo: trim to len_max - end function tList_asStrings From 9119254210a9e186d8eb132627196978980ccbb6 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 17:26:42 +0200 Subject: [PATCH 08/36] materialpoint => homogenization especially as prefix for global data (clear name spaces) --- src/CPFEM.f90 | 14 +++---- src/grid/grid_mech_FEM.f90 | 30 +++++++-------- src/grid/grid_mech_spectral_basic.f90 | 4 +- src/grid/grid_mech_spectral_polarisation.f90 | 6 +-- src/grid/spectral_utilities.f90 | 20 +++++----- src/homogenization.f90 | 40 ++++++++++---------- src/material.f90 | 2 +- src/mesh/FEM_utilities.f90 | 2 +- src/mesh/mesh_mech_FEM.f90 | 20 +++++----- src/results.f90 | 24 ++++++------ 10 files changed, 81 insertions(+), 81 deletions(-) diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index fa774de66..a55b9e4ce 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -184,8 +184,8 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyS temperature(material_homogenizationAt(elCP))%p(thermalMapping(material_homogenizationAt(elCP))%p(ip,elCP)) = & temperature_inp end select chosenThermal1 - materialpoint_F0(1:3,1:3,ip,elCP) = ffn - materialpoint_F(1:3,1:3,ip,elCP) = ffn1 + homogenization_F0(1:3,1:3,ip,elCP) = ffn + homogenization_F(1:3,1:3,ip,elCP) = ffn1 if (iand(mode, CPFEM_CALCRESULTS) /= 0_pInt) then @@ -212,17 +212,17 @@ subroutine CPFEM_general(mode, ffn, ffn1, temperature_inp, dt, elFE, ip, cauchyS else terminalIllness ! translate from P to sigma - Kirchhoff = matmul(materialpoint_P(1:3,1:3,ip,elCP), transpose(materialpoint_F(1:3,1:3,ip,elCP))) - J_inverse = 1.0_pReal / math_det33(materialpoint_F(1:3,1:3,ip,elCP)) + Kirchhoff = matmul(homogenization_P(1:3,1:3,ip,elCP), transpose(homogenization_F(1:3,1:3,ip,elCP))) + J_inverse = 1.0_pReal / math_det33(homogenization_F(1:3,1:3,ip,elCP)) CPFEM_cs(1:6,ip,elCP) = math_sym33to6(J_inverse * Kirchhoff,weighted=.false.) ! translate from dP/dF to dCS/dE H = 0.0_pReal do i=1,3; do j=1,3; do k=1,3; do l=1,3; do m=1,3; do n=1,3 H(i,j,k,l) = H(i,j,k,l) & - + materialpoint_F(j,m,ip,elCP) * materialpoint_F(l,n,ip,elCP) & - * materialpoint_dPdF(i,m,k,n,ip,elCP) & - - math_delta(j,l) * materialpoint_F(i,m,ip,elCP) * materialpoint_P(k,m,ip,elCP) & + + homogenization_F(j,m,ip,elCP) * homogenization_F(l,n,ip,elCP) & + * homogenization_dPdF(i,m,k,n,ip,elCP) & + - math_delta(j,l) * homogenization_F(i,m,ip,elCP) * homogenization_P(k,m,ip,elCP) & + 0.5_pReal * ( Kirchhoff(j,l)*math_delta(i,k) + Kirchhoff(i,k)*math_delta(j,l) & + Kirchhoff(j,k)*math_delta(i,l) + Kirchhoff(i,l)*math_delta(j,k)) enddo; enddo; enddo; enddo; enddo; enddo diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index dfbd3f2f3..bf3d7752d 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -236,7 +236,7 @@ subroutine grid_mech_FEM_init F = spread(spread(spread(math_I3,3,grid(1)),4,grid(2)),5,grid3) endif restartRead - materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent call utilities_updateCoords(F) call utilities_constitutiveResponse(P_current,temp33_Real,C_volAvg,devNull, & ! stress field, stress avg, global average of stiffness and (min+max)/2 F, & ! target F @@ -364,7 +364,7 @@ subroutine grid_mech_FEM_forward(cutBack,guess,timeinc,timeinc_old,loadCaseTime, F_lastInc = F - materialpoint_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) + homogenization_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) endif !-------------------------------------------------------------------------------------------------- @@ -557,9 +557,9 @@ subroutine formResidual(da_local,x_local, & ii = i-xstart+1; jj = j-ystart+1; kk = k-zstart+1 ele = ele + 1 f_elem = matmul(transpose(BMat),transpose(P_current(1:3,1:3,ii,jj,kk)))*detJ + & - matmul(HGMat,x_elem)*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal + matmul(HGMat,x_elem)*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal ctr = 0 do kk = 0, 1; do jj = 0, 1; do ii = 0, 1 ctr = ctr + 1 @@ -636,18 +636,18 @@ subroutine formJacobian(da_local,x_local,Jac_pre,Jac,dummy,ierr) row = col ele = ele + 1 K_ele = 0.0 - K_ele(1 :8 ,1 :8 ) = HGMat*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal - K_ele(9 :16,9 :16) = HGMat*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal - K_ele(17:24,17:24) = HGMat*(materialpoint_dPdF(1,1,1,1,1,ele) + & - materialpoint_dPdF(2,2,2,2,1,ele) + & - materialpoint_dPdF(3,3,3,3,1,ele))/3.0_pReal + K_ele(1 :8 ,1 :8 ) = HGMat*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal + K_ele(9 :16,9 :16) = HGMat*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal + K_ele(17:24,17:24) = HGMat*(homogenization_dPdF(1,1,1,1,1,ele) + & + homogenization_dPdF(2,2,2,2,1,ele) + & + homogenization_dPdF(3,3,3,3,1,ele))/3.0_pReal K_ele = K_ele + & matmul(transpose(BMatFull), & - matmul(reshape(reshape(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,ele), & + matmul(reshape(reshape(homogenization_dPdF(1:3,1:3,1:3,1:3,1,ele), & shape=[3,3,3,3], order=[2,1,4,3]),shape=[9,9]),BMatFull))*detJ call MatSetValuesStencil(Jac,24,row,24,col,K_ele,ADD_VALUES,ierr) CHKERRQ(ierr) diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index d2f6b40da..ec6c3a540 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -198,7 +198,7 @@ subroutine grid_mech_spectral_basic_init F = reshape(F_lastInc,[9,grid(1),grid(2),grid3]) endif restartRead - materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent call utilities_updateCoords(reshape(F,shape(F_lastInc))) call utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F @@ -324,7 +324,7 @@ subroutine grid_mech_spectral_basic_forward(cutBack,guess,timeinc,timeinc_old,lo rotation_BC%rotate(F_aimDot,active=.true.)) F_lastInc = reshape(F,[3,3,grid(1),grid(2),grid3]) - materialpoint_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) + homogenization_F0 = reshape(F, [3,3,1,product(grid(1:2))*grid3]) endif !-------------------------------------------------------------------------------------------------- diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 3d495ddf0..8f9ea81b3 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -224,7 +224,7 @@ subroutine grid_mech_spectral_polarisation_init F_tau_lastInc = 2.0_pReal*F_lastInc endif restartRead - materialpoint_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent + homogenization_F0 = reshape(F_lastInc, [3,3,1,product(grid(1:2))*grid3]) ! set starting condition for materialpoint_stressAndItsTangent call utilities_updateCoords(reshape(F,shape(F_lastInc))) call utilities_constitutiveResponse(P,temp33_Real,C_volAvg,C_minMaxAvg, & ! stress field, stress avg, global average of stiffness and (min+max)/2 reshape(F,shape(F_lastInc)), & ! target F @@ -364,7 +364,7 @@ subroutine grid_mech_spectral_polarisation_forward(cutBack,guess,timeinc,timeinc F_lastInc = reshape(F, [3,3,grid(1),grid(2),grid3]) F_tau_lastInc = reshape(F_tau,[3,3,grid(1),grid(2),grid3]) - materialpoint_F0 = reshape(F,[3,3,1,product(grid(1:2))*grid3]) + homogenization_F0 = reshape(F,[3,3,1,product(grid(1:2))*grid3]) endif !-------------------------------------------------------------------------------------------------- @@ -606,7 +606,7 @@ subroutine formResidual(in, FandF_tau, & do k = 1, grid3; do j = 1, grid(2); do i = 1, grid(1) e = e + 1 residual_F(1:3,1:3,i,j,k) = & - math_mul3333xx33(math_invSym3333(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,e) + C_scale), & + math_mul3333xx33(math_invSym3333(homogenization_dPdF(1:3,1:3,1:3,1:3,1,e) + C_scale), & residual_F(1:3,1:3,i,j,k) - matmul(F(1:3,1:3,i,j,k), & math_mul3333xx33(C_scale,F_tau(1:3,1:3,i,j,k) - F(1:3,1:3,i,j,k) - math_I3))) & + residual_F_tau(1:3,1:3,i,j,k) diff --git a/src/grid/spectral_utilities.f90 b/src/grid/spectral_utilities.f90 index 77047a317..fddd1885f 100644 --- a/src/grid/spectral_utilities.f90 +++ b/src/grid/spectral_utilities.f90 @@ -802,7 +802,7 @@ end subroutine utilities_fourierTensorDivergence !-------------------------------------------------------------------------------------------------- -!> @brief calculate constitutive response from materialpoint_F0 to F during timeinc +!> @brief calculate constitutive response from homogenization_F0 to F during timeinc !-------------------------------------------------------------------------------------------------- subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& F,timeinc,rotation_BC) @@ -824,11 +824,11 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& print'(/,a)', ' ... evaluating constitutive response ......................................' flush(IO_STDOUT) - materialpoint_F = reshape(F,[3,3,1,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field + homogenization_F = reshape(F,[3,3,1,product(grid(1:2))*grid3]) ! set materialpoint target F to estimated field call materialpoint_stressAndItsTangent(timeinc) ! calculate P field - P = reshape(materialpoint_P, [3,3,grid(1),grid(2),grid3]) + P = reshape(homogenization_P, [3,3,grid(1),grid(2),grid3]) P_av = sum(sum(sum(P,dim=5),dim=4),dim=3) * wgt ! average of P call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) if (debugRotation) & @@ -845,13 +845,13 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& dPdF_min = huge(1.0_pReal) dPdF_norm_min = huge(1.0_pReal) do i = 1, product(grid(1:2))*grid3 - if (dPdF_norm_max < sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then - dPdF_max = materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i) - dPdF_norm_max = sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) + if (dPdF_norm_max < sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then + dPdF_max = homogenization_dPdF(1:3,1:3,1:3,1:3,1,i) + dPdF_norm_max = sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) endif - if (dPdF_norm_min > sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then - dPdF_min = materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i) - dPdF_norm_min = sum(materialpoint_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) + if (dPdF_norm_min > sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal)) then + dPdF_min = homogenization_dPdF(1:3,1:3,1:3,1:3,1,i) + dPdF_norm_min = sum(homogenization_dPdF(1:3,1:3,1:3,1:3,1,i)**2.0_pReal) endif end do @@ -869,7 +869,7 @@ subroutine utilities_constitutiveResponse(P,P_av,C_volAvg,C_minmaxAvg,& C_minmaxAvg = 0.5_pReal*(dPdF_max + dPdF_min) - C_volAvg = sum(sum(materialpoint_dPdF,dim=6),dim=5) + C_volAvg = sum(sum(homogenization_dPdF,dim=6),dim=5) call MPI_Allreduce(MPI_IN_PLACE,C_volAvg,81,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) C_volAvg = C_volAvg * wgt diff --git a/src/homogenization.f90 b/src/homogenization.f90 index 75f93dc08..c717c24fd 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -31,12 +31,12 @@ module homogenization !-------------------------------------------------------------------------------------------------- ! General variables for the homogenization at a material point 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 + homogenization_F0, & !< def grad of IP at start of FE increment + homogenization_F !< def grad of IP to be reached at end of FE increment real(pReal), dimension(:,:,:,:), allocatable, public, protected :: & - materialpoint_P !< first P--K stress of IP + homogenization_P !< first P--K stress of IP real(pReal), dimension(:,:,:,:,:,:), allocatable, public, protected :: & - materialpoint_dPdF !< tangent of first P--K stress at IP + homogenization_dPdF !< tangent of first P--K stress at IP type :: tNumerics integer :: & @@ -181,10 +181,10 @@ subroutine homogenization_init !-------------------------------------------------------------------------------------------------- ! allocate and initialize global variables - allocate(materialpoint_dPdF(3,3,3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) - materialpoint_F0 = spread(spread(math_I3,3,discretization_nIP),4,discretization_nElem) ! initialize to identity - materialpoint_F = materialpoint_F0 ! initialize to identity - allocate(materialpoint_P(3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) + allocate(homogenization_dPdF(3,3,3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) + homogenization_F0 = spread(spread(math_I3,3,discretization_nIP),4,discretization_nElem) ! initialize to identity + homogenization_F = homogenization_F0 ! initialize to identity + allocate(homogenization_P(3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) print'(/,a)', ' <<<+- homogenization init -+>>>'; flush(IO_STDOUT) @@ -330,8 +330,8 @@ subroutine materialpoint_stressAndItsTangent(dt) myNgrains = homogenization_Nconstituent(material_homogenizationAt(e)) IpLooping2: do i = FEsolving_execIP(1),FEsolving_execIP(2) if(requested(i,e) .and. .not. doneAndHappy(1,i,e)) then ! requested but not yet done - call partitionDeformation(materialpoint_F0(1:3,1:3,i,e) & - + (materialpoint_F(1:3,1:3,i,e)-materialpoint_F0(1:3,1:3,i,e))& + call partitionDeformation(homogenization_F0(1:3,1:3,i,e) & + + (homogenization_F(1:3,1:3,i,e)-homogenization_F0(1:3,1:3,i,e))& *(subStep(i,e)+subFrac(i,e)), & i,e) crystallite_dt(1:myNgrains,i,e) = dt*subStep(i,e) ! propagate materialpoint dt to grains @@ -357,8 +357,8 @@ subroutine materialpoint_stressAndItsTangent(dt) doneAndHappy(1:2,i,e) = [.true.,.false.] else doneAndHappy(1:2,i,e) = updateState(dt*subStep(i,e), & - materialpoint_F0(1:3,1:3,i,e) & - + (materialpoint_F(1:3,1:3,i,e)-materialpoint_F0(1:3,1:3,i,e)) & + homogenization_F0(1:3,1:3,i,e) & + + (homogenization_F(1:3,1:3,i,e)-homogenization_F0(1:3,1:3,i,e)) & *(subStep(i,e)+subFrac(i,e)), & i,e) converged(i,e) = all(doneAndHappy(1:2,i,e)) ! converged if done and happy @@ -492,16 +492,16 @@ subroutine averageStressAndItsTangent(ip,el) chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) case (HOMOGENIZATION_NONE_ID) chosenHomogenization - materialpoint_P(1:3,1:3,ip,el) = crystallite_P(1:3,1:3,1,ip,el) - materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_stressTangent(1,ip,el) + homogenization_P(1:3,1:3,ip,el) = crystallite_P(1:3,1:3,1,ip,el) + homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_stressTangent(1,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization do c = 1, homogenization_Nconstituent(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_isostrain_averageStressAndItsTangent(& - materialpoint_P(1:3,1:3,ip,el), & - materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& + homogenization_P(1:3,1:3,ip,el), & + homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el),& crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) @@ -511,8 +511,8 @@ subroutine averageStressAndItsTangent(ip,el) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_RGC_averageStressAndItsTangent(& - materialpoint_P(1:3,1:3,ip,el), & - materialpoint_dPdF(1:3,1:3,1:3,1:3,ip,el),& + homogenization_P(1:3,1:3,ip,el), & + homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el),& crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) @@ -539,10 +539,10 @@ subroutine homogenization_results group = trim(group_base)//'/generic' call results_closeGroup(results_addGroup(group)) - !temp = reshape(materialpoint_F,[3,3,discretization_nIP*discretization_nElem]) + !temp = reshape(homogenization_F,[3,3,discretization_nIP*discretization_nElem]) !call results_writeDataset(group,temp,'F',& ! 'deformation gradient','1') - !temp = reshape(materialpoint_P,[3,3,discretization_nIP*discretization_nElem]) + !temp = reshape(homogenization_P,[3,3,discretization_nIP*discretization_nElem]) !call results_writeDataset(group,temp,'P',& ! '1st Piola-Kirchhoff stress','Pa') diff --git a/src/material.f90 b/src/material.f90 index 7197506f4..a9459e42d 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -177,7 +177,7 @@ subroutine material_init(restart) if (.not. restart) then call results_openJobFile call results_mapping_constituent(material_phaseAt,material_phaseMemberAt,material_name_phase) - call results_mapping_materialpoint(material_homogenizationAt,material_homogenizationMemberAt,material_name_homogenization) + call results_mapping_homogenization(material_homogenizationAt,material_homogenizationMemberAt,material_name_homogenization) call results_closeJobFile endif diff --git a/src/mesh/FEM_utilities.f90 b/src/mesh/FEM_utilities.f90 index 4927d0c1c..4c958ee2e 100644 --- a/src/mesh/FEM_utilities.f90 +++ b/src/mesh/FEM_utilities.f90 @@ -164,7 +164,7 @@ subroutine utilities_constitutiveResponse(timeinc,P_av,forwardData) cutBack = .false. ! reset cutBack status - P_av = sum(sum(materialpoint_P,dim=4),dim=3) * wgt ! average of P + P_av = sum(sum(homogenization_P,dim=4),dim=3) * wgt ! average of P call MPI_Allreduce(MPI_IN_PLACE,P_av,9,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD,ierr) end subroutine utilities_constitutiveResponse diff --git a/src/mesh/mesh_mech_FEM.f90 b/src/mesh/mesh_mech_FEM.f90 index f7d33adcf..8aa084ac8 100644 --- a/src/mesh/mesh_mech_FEM.f90 +++ b/src/mesh/mesh_mech_FEM.f90 @@ -400,15 +400,15 @@ subroutine FEM_mech_formResidual(dm_local,xx_local,f_local,dummy,ierr) (((qPt*nBasis + basis)*dimPlex + comp)*dimPlex+comp+1)*dimPlex)) enddo enddo - materialpoint_F(1:dimPlex,1:dimPlex,qPt+1,cell+1) = & + homogenization_F(1:dimPlex,1:dimPlex,qPt+1,cell+1) = & reshape(matmul(BMat,x_scal),shape=[dimPlex,dimPlex], order=[2,1]) enddo if (num%BBarStabilisation) then - detFAvg = math_det33(sum(materialpoint_F(1:3,1:3,1:nQuadrature,cell+1),dim=3)/real(nQuadrature)) + detFAvg = math_det33(sum(homogenization_F(1:3,1:3,1:nQuadrature,cell+1),dim=3)/real(nQuadrature)) do qPt = 1, nQuadrature - materialpoint_F(1:dimPlex,1:dimPlex,qPt,cell+1) = & - materialpoint_F(1:dimPlex,1:dimPlex,qPt,cell+1)* & - (detFAvg/math_det33(materialpoint_F(1:3,1:3,qPt,cell+1)))**(1.0/real(dimPlex)) + homogenization_F(1:dimPlex,1:dimPlex,qPt,cell+1) = & + homogenization_F(1:dimPlex,1:dimPlex,qPt,cell+1)* & + (detFAvg/math_det33(homogenization_F(1:3,1:3,qPt,cell+1)))**(1.0/real(dimPlex)) enddo endif @@ -443,7 +443,7 @@ subroutine FEM_mech_formResidual(dm_local,xx_local,f_local,dummy,ierr) enddo f_scal = f_scal + & matmul(transpose(BMat), & - reshape(transpose(materialpoint_P(1:dimPlex,1:dimPlex,qPt+1,cell+1)), & + reshape(transpose(homogenization_P(1:dimPlex,1:dimPlex,qPt+1,cell+1)), & shape=[dimPlex*dimPlex]))*qWeights(qPt+1) enddo f_scal = f_scal*abs(detJ) @@ -545,7 +545,7 @@ subroutine FEM_mech_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) (((qPt*nBasis + basis)*dimPlex + comp)*dimPlex+comp+1)*dimPlex)) enddo enddo - MatA = matmul(reshape(reshape(materialpoint_dPdF(1:dimPlex,1:dimPlex,1:dimPlex,1:dimPlex,qPt+1,cell+1), & + MatA = matmul(reshape(reshape(homogenization_dPdF(1:dimPlex,1:dimPlex,1:dimPlex,1:dimPlex,qPt+1,cell+1), & shape=[dimPlex,dimPlex,dimPlex,dimPlex], order=[2,1,4,3]), & shape=[dimPlex*dimPlex,dimPlex*dimPlex]),BMat)*qWeights(qPt+1) if (num%BBarStabilisation) then @@ -553,12 +553,12 @@ subroutine FEM_mech_formJacobian(dm_local,xx_local,Jac_pre,Jac,dummy,ierr) FInv = math_inv33(F) K_eA = K_eA + matmul(transpose(BMat),MatA)*math_det33(FInv)**(1.0/real(dimPlex)) K_eB = K_eB - & - matmul(transpose(matmul(reshape(materialpoint_F(1:dimPlex,1:dimPlex,qPt+1,cell+1), & + matmul(transpose(matmul(reshape(homogenization_F(1:dimPlex,1:dimPlex,qPt+1,cell+1), & shape=[dimPlex*dimPlex,1]), & matmul(reshape(FInv(1:dimPlex,1:dimPlex), & shape=[1,dimPlex*dimPlex],order=[2,1]),BMat))),MatA) MatB = MatB + & - matmul(reshape(materialpoint_F(1:dimPlex,1:dimPlex,qPt+1,cell+1),shape=[1,dimPlex*dimPlex]),MatA) + matmul(reshape(homogenization_F(1:dimPlex,1:dimPlex,qPt+1,cell+1),shape=[1,dimPlex*dimPlex]),MatA) FAvg = FAvg + F BMatAvg = BMatAvg + BMat else @@ -630,7 +630,7 @@ subroutine FEM_mech_forward(guess,timeinc,timeinc_old,fieldBC) ! forward last inc if (guess .and. .not. cutBack) then ForwardData = .True. - materialpoint_F0 = materialpoint_F + homogenization_F0 = homogenization_F call SNESGetDM(mech_snes,dm_local,ierr); CHKERRQ(ierr) !< retrieve mesh info from mech_snes into dm_local call DMGetSection(dm_local,section,ierr); CHKERRQ(ierr) call DMGetLocalVector(dm_local,x_local,ierr); CHKERRQ(ierr) diff --git a/src/results.f90 b/src/results.f90 index aec90d7be..1ccc6bfab 100644 --- a/src/results.f90 +++ b/src/results.f90 @@ -56,7 +56,7 @@ module results results_addAttribute, & results_removeLink, & results_mapping_constituent, & - results_mapping_materialpoint + results_mapping_homogenization contains subroutine results_init(restart) @@ -644,7 +644,7 @@ end subroutine results_mapping_constituent !-------------------------------------------------------------------------------------------------- !> @brief adds the unique mapping from spatial position and constituent ID to results !-------------------------------------------------------------------------------------------------- -subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) +subroutine results_mapping_homogenization(homogenizationAt,memberAtLocal,label) integer, dimension(:), intent(in) :: homogenizationAt !< homogenization section at (element) integer, dimension(:,:), intent(in) :: memberAtLocal !< homogenization member at (IP,element) @@ -711,13 +711,13 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) ! MPI settings and communication #ifdef PETSc call h5pset_dxpl_mpio_f(plist_id, H5FD_MPIO_COLLECTIVE_F, ierr) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5pset_dxpl_mpio_f') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5pset_dxpl_mpio_f') call MPI_allreduce(MPI_IN_PLACE,writeSize,worldsize,MPI_INT,MPI_SUM,PETSC_COMM_WORLD,ierr) ! get output at each process - if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_materialpoint: MPI_allreduce/writeSize') + if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_homogenization: MPI_allreduce/writeSize') call MPI_allreduce(MPI_IN_PLACE,memberOffset,size(memberOffset),MPI_INT,MPI_SUM,PETSC_COMM_WORLD,ierr)! get offset at each process - if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_materialpoint: MPI_allreduce/memberOffset') + if (ierr /= 0) call IO_error(894,ext_msg='results_mapping_homogenization: MPI_allreduce/memberOffset') #endif myShape = int([writeSize(worldrank)], HSIZE_T) @@ -727,13 +727,13 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) !-------------------------------------------------------------------------------------------------- ! create dataspace in memory (local shape = hyperslab) and in file (global shape) call h5screate_simple_f(1,myShape,memspace_id,ierr,myShape) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5screate_simple_f/memspace_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5screate_simple_f/memspace_id') call h5screate_simple_f(1,totalShape,filespace_id,ierr,totalShape) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5screate_simple_f/filespace_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5screate_simple_f/filespace_id') call h5sselect_hyperslab_f(filespace_id, H5S_SELECT_SET_F, myOffset, myShape, ierr) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5sselect_hyperslab_f') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5sselect_hyperslab_f') !--------------------------------------------------------------------------------------------------- ! expand phaseAt to consider IPs (is not stored per IP) @@ -753,14 +753,14 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) loc_id = results_openGroup('/mapping') call h5dcreate_f(loc_id, 'homogenization', dtype_id, filespace_id, dset_id, ierr) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5dcreate_f') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5dcreate_f') call h5dwrite_f(dset_id, name_id, reshape(label(pack(homogenizationAtMaterialpoint,.true.)),myShape), & myShape, ierr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5dwrite_f/name_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5dwrite_f/name_id') call h5dwrite_f(dset_id, position_id, reshape(pack(memberAtGlobal,.true.),myShape), & myShape, ierr, file_space_id = filespace_id, mem_space_id = memspace_id, xfer_prp = plist_id) - if (ierr < 0) call IO_error(1,ext_msg='results_mapping_materialpoint: h5dwrite_f/position_id') + if (ierr < 0) call IO_error(1,ext_msg='results_mapping_homogenization: h5dwrite_f/position_id') !-------------------------------------------------------------------------------------------------- ! close all @@ -776,7 +776,7 @@ subroutine results_mapping_materialpoint(homogenizationAt,memberAtLocal,label) ! for backward compatibility call results_setLink('/mapping/homogenization','/mapping/cellResults/materialpoint') -end subroutine results_mapping_materialpoint +end subroutine results_mapping_homogenization !-------------------------------------------------------------------------------------------------- From fd4cdf965b189e3bbb4b76f7c85067fd5111362e Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 24 Oct 2020 18:00:17 +0200 Subject: [PATCH 09/36] modularized --- src/material.f90 | 85 +++++++++++++++++++++++++++++------------------- 1 file changed, 52 insertions(+), 33 deletions(-) diff --git a/src/material.f90 b/src/material.f90 index a9459e42d..58a5c4080 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -149,7 +149,7 @@ contains subroutine material_init(restart) logical, intent(in) :: restart - integer :: myHomog + integer :: myHomog print'(/,a)', ' <<<+- material init -+>>>'; flush(IO_STDOUT) @@ -298,48 +298,22 @@ subroutine material_parseMaterial counterPhase, & counterHomogenization - character(len=pStringLen) :: sectionName - real(pReal) :: & frac integer :: & e, i, c, & - p, h, & - m + h materials => config_material%get('material') phases => config_material%get('phase') homogenizations => config_material%get('homogenization') + call sanityCheck(materials, homogenizations) + material_name_phase = getKeys(phases) + material_name_homogenization = getKeys(homogenizations) - if(any(discretization_materialAt > materials%length)) & - call IO_error(155,ext_msg='More materials requested than found in material.yaml') - -!-------------------------------------------------------------------------------------------------- -! determine name of phases/homogenizations - allocate(material_name_phase(phases%length)) - do p = 1, phases%length - write(sectionName,'(i0,a)') p,'_' - material_name_phase(p) = trim(adjustl(sectionName))//phases%getKey(p) !ToDo: remove prefix - enddo - - allocate(material_name_homogenization(homogenizations%length)) - do h = 1, homogenizations%length - write(sectionName,'(i0,a)') h,'_' - material_name_homogenization(h) = trim(adjustl(sectionName))//homogenizations%getKey(h) !ToDo: remove prefix - enddo - -!-------------------------------------------------------------------------------------------------- -! sanity check: Matching # of constituents - do m = 1, materials%length - material => materials%get(m) - constituents => material%get('constituents') - homogenization => homogenizations%get(material%get_asString('homogenization')) - if(constituents%length /= homogenization%get_asInt('N_constituents')) call IO_error(148) - enddo - - allocate(homogenization_Nconstituent(size(material_name_homogenization))) - do h=1, size(material_name_homogenization) + allocate(homogenization_Nconstituent(homogenizations%length)) + do h=1, homogenizations%length homogenization => homogenizations%get(h) homogenization_Nconstituent(h) = homogenization%get_asInt('N_constituents') enddo @@ -386,4 +360,49 @@ subroutine material_parseMaterial end subroutine material_parseMaterial +!-------------------------------------------------------------------------------------------------- +!> @brief Check if material.yaml is consistent and contains sufficient # of materials +!-------------------------------------------------------------------------------------------------- +subroutine sanityCheck(materials,homogenizations) + + class(tNode), intent(in) :: materials, & + homogenizations + + class(tNode), pointer :: material, & + homogenization, & + constituents + integer :: m + + if(maxval(discretization_materialAt) > materials%length) & + call IO_error(155,ext_msg='More materials requested than found in material.yaml') + + do m = 1, materials%length + material => materials%get(m) + constituents => material%get('constituents') + homogenization => homogenizations%get(material%get_asString('homogenization')) + if(constituents%length /= homogenization%get_asInt('N_constituents')) call IO_error(148) + enddo + +end subroutine sanityCheck + + +!-------------------------------------------------------------------------------------------------- +!> @brief Get all keys from a dictionary (currently with #_ prefix) +!-------------------------------------------------------------------------------------------------- +function getKeys(dict) + + class(tNode), intent(in) :: dict + character(len=pStringLen), dimension(:), allocatable :: getKeys + + integer :: i + character(len=pStringLen) :: sectionName + + allocate(getKeys(dict%length)) + do i=1, dict%length + write(sectionName,'(i0,a)') i,'_' + getKeys(i) = trim(adjustl(sectionName))//dict%getKey(i) !ToDo: remove prefix + enddo + +end function getKeys + end module material From 93a791abd19ff1f8143633f06bba073d088fbda0 Mon Sep 17 00:00:00 2001 From: Test User Date: Sat, 24 Oct 2020 19:17:53 +0200 Subject: [PATCH 10/36] [skip ci] updated version information after successful test of v3.0.0-alpha-556-gf379aecd8 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 64636c48f..2e3222238 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha-553-g5eee0d74f +v3.0.0-alpha-556-gf379aecd8 From 5834d9501231f877ada4c065b5e43ff1aa68cae7 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 26 Oct 2020 21:38:24 +0100 Subject: [PATCH 11/36] support for more types - allow to directly use Config and its sublasses (cast to dict) - convert numpy arrays --- python/damask/_config.py | 15 ++++++++++++++- 1 file changed, 14 insertions(+), 1 deletion(-) diff --git a/python/damask/_config.py b/python/damask/_config.py index c9130d7aa..40b9e918f 100644 --- a/python/damask/_config.py +++ b/python/damask/_config.py @@ -1,6 +1,7 @@ from io import StringIO import abc +import numpy as np import yaml class NiceDumper(yaml.SafeDumper): @@ -15,6 +16,11 @@ class NiceDumper(yaml.SafeDumper): def increase_indent(self, flow=False, indentless=False): return super().increase_indent(flow, False) + def represent_data(self, data): + """Cast Config objects and its qsubclasses to dict.""" + return self.represent_data(dict(data)) if isinstance(data, dict) and type(data) != dict else \ + super().represent_data(data) + class Config(dict): """YAML-based configuration.""" @@ -64,7 +70,14 @@ class Config(dict): kwargs['width'] = 256 if 'default_flow_style' not in kwargs: kwargs['default_flow_style'] = None - fhandle.write(yaml.dump(dict(self),Dumper=NiceDumper,**kwargs)) + + def array_representer(dumper, data): + """Convert numpy array to list of native types.""" + return dumper.represent_list([d.item() for d in data]) + + NiceDumper.add_representer(np.ndarray, array_representer) + + fhandle.write(yaml.dump(self,Dumper=NiceDumper,**kwargs)) @property From 4159ae5f3d6ffad63fb71b98c79d385953681523 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Mon, 26 Oct 2020 21:43:21 +0100 Subject: [PATCH 12/36] keep order of keys --- python/damask/_config.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/python/damask/_config.py b/python/damask/_config.py index 40b9e918f..f2ffbbbf3 100644 --- a/python/damask/_config.py +++ b/python/damask/_config.py @@ -17,7 +17,7 @@ class NiceDumper(yaml.SafeDumper): return super().increase_indent(flow, False) def represent_data(self, data): - """Cast Config objects and its qsubclasses to dict.""" + """Cast Config objects and its subclasses to dict.""" return self.represent_data(dict(data)) if isinstance(data, dict) and type(data) != dict else \ super().represent_data(data) @@ -70,6 +70,8 @@ class Config(dict): kwargs['width'] = 256 if 'default_flow_style' not in kwargs: kwargs['default_flow_style'] = None + if 'sort_keys' not in kwargs: + kwargs['sort_keys'] = False def array_representer(dumper, data): """Convert numpy array to list of native types.""" From 4dd387d7f51c5b25535126f2c8733a3a91f25bec Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 06:34:05 +0100 Subject: [PATCH 13/36] compatibility with old pyyaml + test for numpy --- python/damask/_config.py | 6 +++++- python/tests/test_Config.py | 3 +++ 2 files changed, 8 insertions(+), 1 deletion(-) diff --git a/python/damask/_config.py b/python/damask/_config.py index f2ffbbbf3..24245f4bd 100644 --- a/python/damask/_config.py +++ b/python/damask/_config.py @@ -79,7 +79,11 @@ class Config(dict): NiceDumper.add_representer(np.ndarray, array_representer) - fhandle.write(yaml.dump(self,Dumper=NiceDumper,**kwargs)) + try: + fhandle.write(yaml.dump(self,Dumper=NiceDumper,**kwargs)) + except TypeError: # compatibility with old pyyaml + del kwargs['sort_keys'] + fhandle.write(yaml.dump(self,Dumper=NiceDumper,**kwargs)) @property diff --git a/python/tests/test_Config.py b/python/tests/test_Config.py index e715ad763..67c419b3e 100644 --- a/python/tests/test_Config.py +++ b/python/tests/test_Config.py @@ -1,4 +1,5 @@ import pytest +import numpy as np from damask import Config @@ -29,6 +30,8 @@ class TestConfig: f.write(config.__repr__()) assert Config.load(tmp_path/'config.yaml') == config + def test_numpy(self,tmp_path): + assert Config({'A':np.ones(3,'i')}).__repr__() == Config({'A':[1,1,1]}).__repr__() def test_abstract_is_valid(self): assert Config().is_valid is None From 201a62d7c9c882f10a0c15f393f613b56ecc7d89 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 13:19:53 +0100 Subject: [PATCH 14/36] easier visualization if verts have a point-cell --- python/damask/_vtk.py | 16 ++++++++++++---- python/tests/reference/VTK/polyData.vtp | 10 +++++----- 2 files changed, 17 insertions(+), 9 deletions(-) diff --git a/python/damask/_vtk.py b/python/damask/_vtk.py index c1fe52f38..546ff4c88 100644 --- a/python/damask/_vtk.py +++ b/python/damask/_vtk.py @@ -108,11 +108,18 @@ class VTK: Spatial position of the points. """ + N = points.shape[0] vtk_points = vtk.vtkPoints() vtk_points.SetData(np_to_vtk(points)) + vtk_cells = vtk.vtkCellArray() + vtk_cells.SetNumberOfCells(N) + vtk_cells.SetCells(N,np_to_vtkIdTypeArray(np.stack((np.ones (N,dtype=np.int64), + np.arange(N,dtype=np.int64)),axis=1).ravel(),deep=True)) + vtk_data = vtk.vtkPolyData() vtk_data.SetPoints(vtk_points) + vtk_data.SetVerts(vtk_cells) return VTK(vtk_data) @@ -164,6 +171,7 @@ class VTK: return VTK(vtk_data) + @staticmethod def _write(writer): """Wrapper for parallel writing.""" @@ -192,7 +200,7 @@ class VTK: default_ext = writer.GetDefaultFileExtension() ext = Path(fname).suffix if ext and ext != '.'+default_ext: - raise ValueError(f'Given extension {ext} does not match default .{default_ext}') + raise ValueError(f'Given extension "{ext}" does not match default ".{default_ext}"') writer.SetFileName(str(Path(fname).with_suffix('.'+default_ext))) if compress: writer.SetCompressorTypeToZLib() @@ -238,10 +246,10 @@ class VTK: else data).reshape(N_data,-1),deep=True) # avoid large files d.SetName(label) - if N_data == N_cells: - self.vtk_data.GetCellData().AddArray(d) - elif N_data == N_points: + if N_data == N_points: self.vtk_data.GetPointData().AddArray(d) + elif N_data == N_cells: + self.vtk_data.GetCellData().AddArray(d) else: raise ValueError(f'Cell / point count ({N_cells} / {N_points}) differs from data ({N_data}).') elif isinstance(data,pd.DataFrame): diff --git a/python/tests/reference/VTK/polyData.vtp b/python/tests/reference/VTK/polyData.vtp index 6ed05f67f..dc4b5f149 100644 --- a/python/tests/reference/VTK/polyData.vtp +++ b/python/tests/reference/VTK/polyData.vtp @@ -1,7 +1,7 @@ - + AQAAAACAAAB4AAAAVgAAAA==eF5jYICBhv2WfY9tLfuS7Ypk3PeDaCDf7okF3/7Vq1bZrV6lZQ+k94HEgHL2QHovUM7+iUUfiG0LlQdhkH77Ipnj9iB5qFp7kBjQDiBmcADRANsaLXM= @@ -31,11 +31,11 @@ - - AAAAAACAAAAAAAAA + + AQAAAACAAABQAAAAIgAAAA==eF4txbcBACAIADAsiP7/sAPJkog2PL28nT4uXz9/BXgALg== - - AAAAAACAAAAAAAAA + + AQAAAACAAABQAAAAIgAAAA==eF4txbcBACAIADA76v8HM5As6a0MTy9vH4evn78TBzAAOA== From fa67a2ddf8d5ec8e555e4d3a1997691d6c4d36a2 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 13:35:14 +0100 Subject: [PATCH 15/36] cell_coordinates is now a property --- python/damask/_result.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/python/damask/_result.py b/python/damask/_result.py index 5e8a9a9d0..ab0c5f13e 100644 --- a/python/damask/_result.py +++ b/python/damask/_result.py @@ -1219,7 +1219,7 @@ class Result: f['/geometry/T_c'][()]-1, f['/geometry/T_c'].attrs['VTK_TYPE'].decode()) elif mode.lower()=='point': - v = VTK.from_polyData(self.cell_coordinates()) + v = VTK.from_polyData(self.cell_coordinates) N_digits = int(np.floor(np.log10(max(1,int(self.increments[-1][3:])))))+1 From 0804e6ee74dbb349a2791d6e466c1c640623c4bd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 13:42:49 +0100 Subject: [PATCH 16/36] pythonic names --- python/damask/_geom.py | 4 ++-- python/damask/_result.py | 10 +++++----- python/damask/_vtk.py | 6 +++--- python/tests/test_Geom.py | 2 +- python/tests/test_VTK.py | 14 +++++++------- 5 files changed, 18 insertions(+), 18 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 6366f06ab..f7ac6c437 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -473,7 +473,7 @@ class Geom: Compress with zlib algorithm. Defaults to True. """ - v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) + v = VTK.from_rectilinear_grid(self.grid,self.size,self.origin) v.add(self.material.flatten(order='F'),'material') v.add_comments(self.comments) @@ -508,7 +508,7 @@ class Geom: def show(self): """Show on screen.""" - v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) + v = VTK.from_rectilinear_grid(self.grid,self.size,self.origin) v.show() diff --git a/python/damask/_result.py b/python/damask/_result.py index ab0c5f13e..3f4ca735c 100644 --- a/python/damask/_result.py +++ b/python/damask/_result.py @@ -1212,14 +1212,14 @@ class Result: if mode.lower()=='cell': if self.structured: - v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) + v = VTK.from_rectilinear_grid(self.grid,self.size,self.origin) else: with h5py.File(self.fname,'r') as f: - v = VTK.from_unstructuredGrid(f['/geometry/x_n'][()], - f['/geometry/T_c'][()]-1, - f['/geometry/T_c'].attrs['VTK_TYPE'].decode()) + v = VTK.from_unstructured_grid(f['/geometry/x_n'][()], + f['/geometry/T_c'][()]-1, + f['/geometry/T_c'].attrs['VTK_TYPE'].decode()) elif mode.lower()=='point': - v = VTK.from_polyData(self.cell_coordinates) + v = VTK.from_poly_data(self.cell_coordinates) N_digits = int(np.floor(np.log10(max(1,int(self.increments[-1][3:])))))+1 diff --git a/python/damask/_vtk.py b/python/damask/_vtk.py index 546ff4c88..10f3b482b 100644 --- a/python/damask/_vtk.py +++ b/python/damask/_vtk.py @@ -36,7 +36,7 @@ class VTK: @staticmethod - def from_rectilinearGrid(grid,size,origin=np.zeros(3)): + def from_rectilinear_grid(grid,size,origin=np.zeros(3)): """ Create VTK of type vtk.vtkRectilinearGrid. @@ -64,7 +64,7 @@ class VTK: @staticmethod - def from_unstructuredGrid(nodes,connectivity,cell_type): + def from_unstructured_grid(nodes,connectivity,cell_type): """ Create VTK of type vtk.vtkUnstructuredGrid. @@ -96,7 +96,7 @@ class VTK: @staticmethod - def from_polyData(points): + def from_poly_data(points): """ Create VTK of type vtk.polyData. diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 0f4c1a952..368d8ec8f 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -49,7 +49,7 @@ class TestGeom: assert geom_equal(new,default) def test_invalid_vtr(self,tmp_path): - v = VTK.from_rectilinearGrid(np.random.randint(5,10,3)*2,np.random.random(3) + 1.0) + v = VTK.from_rectilinear_grid(np.random.randint(5,10,3)*2,np.random.random(3) + 1.0) v.save(tmp_path/'no_materialpoint.vtr') with pytest.raises(ValueError): Geom.load(tmp_path/'no_materialpoint.vtr') diff --git a/python/tests/test_VTK.py b/python/tests/test_VTK.py index 81c9eb772..63adb6454 100644 --- a/python/tests/test_VTK.py +++ b/python/tests/test_VTK.py @@ -18,7 +18,7 @@ def default(): """Simple VTK.""" grid = np.array([5,6,7],int) size = np.array([.6,1.,.5]) - return VTK.from_rectilinearGrid(grid,size) + return VTK.from_rectilinear_grid(grid,size) class TestVTK: @@ -30,7 +30,7 @@ class TestVTK: grid = np.random.randint(5,10,3)*2 size = np.random.random(3) + 1.0 origin = np.random.random(3) - v = VTK.from_rectilinearGrid(grid,size,origin) + v = VTK.from_rectilinear_grid(grid,size,origin) string = v.__repr__() v.save(tmp_path/'rectilinearGrid',False) vtr = VTK.load(tmp_path/'rectilinearGrid.vtr') @@ -41,7 +41,7 @@ class TestVTK: def test_polyData(self,tmp_path): points = np.random.rand(100,3) - v = VTK.from_polyData(points) + v = VTK.from_poly_data(points) string = v.__repr__() v.save(tmp_path/'polyData',False) vtp = VTK.load(tmp_path/'polyData.vtp') @@ -60,7 +60,7 @@ class TestVTK: def test_unstructuredGrid(self,tmp_path,cell_type,n): nodes = np.random.rand(n,3) connectivity = np.random.choice(np.arange(n),n,False).reshape(-1,n) - v = VTK.from_unstructuredGrid(nodes,connectivity,cell_type) + v = VTK.from_unstructured_grid(nodes,connectivity,cell_type) string = v.__repr__() v.save(tmp_path/'unstructuredGrid',False) vtu = VTK.load(tmp_path/'unstructuredGrid.vtu') @@ -72,7 +72,7 @@ class TestVTK: def test_parallel_out(self,tmp_path): points = np.random.rand(102,3) - v = VTK.from_polyData(points) + v = VTK.from_poly_data(points) fname_s = tmp_path/'single.vtp' fname_p = tmp_path/'parallel.vtp' v.save(fname_s,False) @@ -121,7 +121,7 @@ class TestVTK: def test_compare_reference_polyData(self,update,reference_dir,tmp_path): points=np.dstack((np.linspace(0.,1.,10),np.linspace(0.,2.,10),np.linspace(-1.,1.,10))).squeeze() - polyData = VTK.from_polyData(points) + polyData = VTK.from_poly_data(points) polyData.add(points,'coordinates') if update: polyData.save(reference_dir/'polyData') @@ -133,7 +133,7 @@ class TestVTK: def test_compare_reference_rectilinearGrid(self,update,reference_dir,tmp_path): grid = np.array([5,6,7],int) size = np.array([.6,1.,.5]) - rectilinearGrid = VTK.from_rectilinearGrid(grid,size) + rectilinearGrid = VTK.from_rectilinear_grid(grid,size) c = grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') n = grid_filters.node_coord0(grid,size).reshape(-1,3,order='F') rectilinearGrid.add(c,'cell') From a5a127b78721ee1a559b278397afbc6bbb026149 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 16:46:08 +0100 Subject: [PATCH 17/36] new tests+ updated results --- .../tests/reference/Result/12grains6x7x8.geom | 61 ------------------- .../tests/reference/Result/12grains6x7x8.vtr | 30 +++++++++ .../Result/6grains6x7x8_single_phase.geom | 61 ------------------- .../Result/6grains6x7x8_single_phase.vtr | 30 +++++++++ python/tests/test_Result.py | 5 ++ 5 files changed, 65 insertions(+), 122 deletions(-) delete mode 100644 python/tests/reference/Result/12grains6x7x8.geom create mode 100644 python/tests/reference/Result/12grains6x7x8.vtr delete mode 100644 python/tests/reference/Result/6grains6x7x8_single_phase.geom create mode 100644 python/tests/reference/Result/6grains6x7x8_single_phase.vtr diff --git a/python/tests/reference/Result/12grains6x7x8.geom b/python/tests/reference/Result/12grains6x7x8.geom deleted file mode 100644 index 4e6836bb8..000000000 --- a/python/tests/reference/Result/12grains6x7x8.geom +++ /dev/null @@ -1,61 +0,0 @@ -4 header -grid a 6 b 7 c 8 -size x 0.75 y 0.875 z 1.0 -origin x 0.0 y 0.0 z 0.0 -homogenization 1 - 9 3 3 10 9 9 - 9 1 1 1 9 9 - 9 11 1 1 7 9 - 7 11 11 7 7 7 - 7 11 11 7 7 7 -12 3 3 10 7 12 -12 3 3 10 10 12 -12 3 3 1 9 9 - 9 1 1 1 9 9 - 9 1 1 1 7 7 - 7 1 1 7 7 7 -12 12 3 7 7 7 -12 3 3 3 12 12 -12 3 3 3 12 12 -12 3 3 1 1 12 - 9 1 1 1 1 9 - 6 1 1 1 8 8 - 7 6 8 8 8 8 -12 12 8 8 8 12 -12 3 3 3 12 12 -12 3 3 3 12 12 - 5 6 6 6 1 12 - 6 6 6 6 8 8 - 6 6 6 8 8 8 - 8 6 8 8 8 8 -12 5 8 8 8 8 -12 5 5 8 8 12 - 5 5 5 3 12 12 - 5 5 6 6 6 5 - 6 6 6 6 6 6 - 6 6 6 6 8 8 - 4 4 6 8 8 8 - 4 4 2 2 2 8 - 5 5 5 2 2 2 - 5 5 5 5 2 5 - 5 5 5 10 10 5 - 6 6 6 6 10 4 - 4 4 11 11 2 4 - 4 4 11 2 2 4 - 4 4 2 2 2 2 - 5 5 5 2 2 2 - 5 5 5 10 10 5 - 5 5 10 10 10 9 - 4 11 11 11 10 9 - 4 4 11 11 11 4 - 4 4 11 11 2 4 - 4 4 2 2 2 2 - 5 5 2 2 2 2 - 5 5 10 10 10 10 - 9 10 10 10 10 9 - 9 11 11 10 9 9 - 4 11 11 11 9 9 - 4 11 11 11 7 7 - 4 4 11 2 7 7 -12 10 10 10 10 7 - 9 10 10 10 10 9 diff --git a/python/tests/reference/Result/12grains6x7x8.vtr b/python/tests/reference/Result/12grains6x7x8.vtr new file mode 100644 index 000000000..02ba9d4e6 --- /dev/null +++ b/python/tests/reference/Result/12grains6x7x8.vtr @@ -0,0 +1,30 @@ + + + + + + AQAAAACAAAARAAAAGQAAAA==eF7LyM/NT0/Ny6xKLMnMz1MwZAAAPsIGPQ== + + + + + + + + AQAAAACAAACACgAA2wAAAA==eF6tlssOgkAQBNcHiPj//2uM25eOZc8oc6kYagzdISzbeM/ZeJ/cgDTk7x/c16zm6fduXAOr/mOS8rqXfDH5mqP6pKHcVY9yJN/ziv5/R/k+8lI/GnnLV2uMm1G5tefXnZ6j6v3bD/nXyQWokU8e5a96tJc8z9H1aY88MfWZek3Xf6XnuBhTr+6fgPKpH9oj3/eS5+/bap/yPafo54buJ/mek3zqJeXu+tRP8vycp15E8tNe6rPaf7fPrk/9eO6q598/1GO1/67v53V6nul8T3n9Oy758p47Sgdl + + + + + AQAAAACAAAA4AAAAHAAAAA==eF5jYEAGB+wh9AUofQNKP4DST6D0C3sAisIGjw== + + + AQAAAACAAABAAAAAHwAAAA==eF5jYEAGB+wh9AUofQNKP4DST6D0Cyj9xh4AwVEHug== + + + AQAAAACAAABIAAAAIgAAAA==eF5jYEAGB+wh9AUofQNKP4DST6D0Cyj9Bkp/sAcAAU8I6Q== + + + + + diff --git a/python/tests/reference/Result/6grains6x7x8_single_phase.geom b/python/tests/reference/Result/6grains6x7x8_single_phase.geom deleted file mode 100644 index 14cae0973..000000000 --- a/python/tests/reference/Result/6grains6x7x8_single_phase.geom +++ /dev/null @@ -1,61 +0,0 @@ -4 header -grid a 6 b 7 c 8 -size x 0.75 y 0.875 z 1.0 -origin x 0.0 y 0.0 z 0.0 -homogenization 1 -3 3 3 4 3 3 -3 1 1 1 3 3 -3 5 1 1 1 3 -1 5 5 1 1 1 -1 5 5 1 1 1 -6 3 3 4 1 6 -6 3 3 4 4 6 -6 3 3 1 3 3 -3 1 1 1 3 3 -3 1 1 1 1 1 -1 1 1 1 1 1 -6 6 3 1 1 1 -6 3 3 3 6 6 -6 3 3 3 6 6 -6 3 3 1 1 6 -3 1 1 1 1 3 -6 1 1 1 2 2 -1 6 2 2 2 2 -6 6 2 2 2 6 -6 3 3 3 6 6 -6 3 3 3 6 6 -5 6 6 6 1 6 -6 6 6 6 2 2 -6 6 6 2 2 2 -2 6 2 2 2 2 -6 5 2 2 2 2 -6 5 5 2 2 6 -5 5 5 3 6 6 -5 5 6 6 6 5 -6 6 6 6 6 6 -6 6 6 6 2 2 -4 4 6 2 2 2 -4 4 2 2 2 2 -5 5 5 2 2 2 -5 5 5 5 2 5 -5 5 5 4 4 5 -6 6 6 6 4 4 -4 4 5 5 2 4 -4 4 5 2 2 4 -4 4 2 2 2 2 -5 5 5 2 2 2 -5 5 5 4 4 5 -5 5 4 4 4 3 -4 5 5 5 4 3 -4 4 5 5 5 4 -4 4 5 5 2 4 -4 4 2 2 2 2 -5 5 2 2 2 2 -5 5 4 4 4 4 -3 4 4 4 4 3 -3 5 5 4 3 3 -4 5 5 5 3 3 -4 5 5 5 1 1 -4 4 5 2 1 1 -6 4 4 4 4 1 -3 4 4 4 4 3 diff --git a/python/tests/reference/Result/6grains6x7x8_single_phase.vtr b/python/tests/reference/Result/6grains6x7x8_single_phase.vtr new file mode 100644 index 000000000..9a7b7dd7d --- /dev/null +++ b/python/tests/reference/Result/6grains6x7x8_single_phase.vtr @@ -0,0 +1,30 @@ + + + + + + AQAAAACAAAARAAAAGQAAAA==eF7LyM/NT0/Ny6xKLMnMz1MwZAAAPsIGPQ== + + + + + + + + AQAAAACAAACACgAAwAAAAA==eF69lcsOwjAQA6GU//9lDsUSGjHyBgq+WGpm09jqY7sc2uA3uR43Gb+/YV/FfXd405S/P93ykmt8vPHRWX3+SpbDZHnj3O8snpqeN+L9TFd4lDmu05ljyn3bj/F5P2yfyNZbnilnc1MuOVZ5mzMu3vpsvbb1T5057Ltk/ZBvfVo/qzznGsdzTvvknO3D8zS+9fjvPlsu4/ifn87bf7DNNf7sPlf59rxYbuPircdp/6s81Z5navoeRav9PACxsANv + + + + + AQAAAACAAAA4AAAAHAAAAA==eF5jYEAGB+wh9AUofQNKP4DST6D0C3sAisIGjw== + + + AQAAAACAAABAAAAAHwAAAA==eF5jYEAGB+wh9AUofQNKP4DST6D0Cyj9xh4AwVEHug== + + + AQAAAACAAABIAAAAIgAAAA==eF5jYEAGB+wh9AUofQNKP4DST6D0Cyj9Bkp/sAcAAU8I6Q== + + + + + diff --git a/python/tests/test_Result.py b/python/tests/test_Result.py index 68b72badf..5ef17cd26 100644 --- a/python/tests/test_Result.py +++ b/python/tests/test_Result.py @@ -341,6 +341,11 @@ class TestResult: os.chdir(tmp_path) default.save_vtk(output) + @pytest.mark.parametrize('mode',['point','cell']) + def test_vtk_mode(self,tmp_path,single_phase,mode): + os.chdir(tmp_path) + single_phase.save_vtk(mode=mode) + def test_XDMF(self,tmp_path,single_phase): os.chdir(tmp_path) single_phase.save_XDMF() From 839be90943bd7fd9edb50f63ae5995cae1a3c092 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 21:27:26 +0100 Subject: [PATCH 18/36] consistent naming Nxxxs => the number of xxx --- src/CPFEM.f90 | 6 +-- src/constitutive.f90 | 6 +-- src/constitutive_damage.f90 | 2 +- src/constitutive_plastic_disloTungsten.f90 | 2 +- src/constitutive_plastic_dislotwin.f90 | 2 +- src/constitutive_plastic_isotropic.f90 | 2 +- src/constitutive_plastic_kinehardening.f90 | 2 +- src/constitutive_plastic_none.f90 | 2 +- src/constitutive_plastic_nonlocal.f90 | 16 ++++---- src/constitutive_plastic_phenopowerlaw.f90 | 2 +- src/constitutive_thermal.f90 | 2 +- src/crystallite.f90 | 38 +++++++++--------- src/damage_local.f90 | 4 +- src/damage_nonlocal.f90 | 12 +++--- src/discretization.f90 | 8 ++-- src/homogenization.f90 | 46 +++++++++++----------- src/homogenization_mech_RGC.f90 | 2 +- src/homogenization_mech_isostrain.f90 | 2 +- src/homogenization_mech_none.f90 | 2 +- src/material.f90 | 28 ++++++------- src/source_damage_anisoBrittle.f90 | 2 +- src/source_damage_anisoDuctile.f90 | 2 +- src/source_damage_isoBrittle.f90 | 2 +- src/source_damage_isoDuctile.f90 | 2 +- src/source_thermal_dissipation.f90 | 2 +- src/source_thermal_externalheat.f90 | 2 +- src/thermal_adiabatic.f90 | 12 +++--- src/thermal_conduction.f90 | 16 ++++---- 28 files changed, 113 insertions(+), 113 deletions(-) diff --git a/src/CPFEM.f90 b/src/CPFEM.f90 index a55b9e4ce..c49ecbcb6 100644 --- a/src/CPFEM.f90 +++ b/src/CPFEM.f90 @@ -107,9 +107,9 @@ subroutine CPFEM_init print'(/,a)', ' <<<+- CPFEM init -+>>>'; flush(IO_STDOUT) - allocate(CPFEM_cs( 6,discretization_nIP,discretization_nElem), source= 0.0_pReal) - allocate(CPFEM_dcsdE( 6,6,discretization_nIP,discretization_nElem), source= 0.0_pReal) - allocate(CPFEM_dcsdE_knownGood(6,6,discretization_nIP,discretization_nElem), source= 0.0_pReal) + allocate(CPFEM_cs( 6,discretization_nIPs,discretization_Nelems), source= 0.0_pReal) + allocate(CPFEM_dcsdE( 6,6,discretization_nIPs,discretization_Nelems), source= 0.0_pReal) + allocate(CPFEM_dcsdE_knownGood(6,6,discretization_nIPs,discretization_Nelems), source= 0.0_pReal) !------------------------------------------------------------------------------ ! read debug options diff --git a/src/constitutive.f90 b/src/constitutive.f90 index 9bbb9d285..c7658b77f 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -127,7 +127,7 @@ module constitutive instance,of,ip,el) real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress - real(pReal), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems), intent(in) :: & F, & !< deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -218,7 +218,7 @@ module constitutive instance, & i, & e - type(rotation), dimension(1,discretization_nIP,discretization_nElem), intent(in) :: & + type(rotation), dimension(1,discretization_nIPs,discretization_Nelems), intent(in) :: & orientation !< crystal orientation end subroutine plastic_nonlocal_updateCompatibility @@ -753,7 +753,7 @@ function constitutive_collectDotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el of real(pReal), intent(in) :: & subdt !< timestep - real(pReal), intent(in), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem) :: & + real(pReal), intent(in), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems) :: & FArray, & !< elastic deformation gradient FpArray !< plastic deformation gradient real(pReal), intent(in), dimension(3,3) :: & diff --git a/src/constitutive_damage.f90 b/src/constitutive_damage.f90 index 39b8534fe..56198987e 100644 --- a/src/constitutive_damage.f90 +++ b/src/constitutive_damage.f90 @@ -184,7 +184,7 @@ module subroutine constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi phiDot = 0.0_pReal dPhiDot_dPhi = 0.0_pReal - do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) phase = material_phaseAt(grain,el) constituent = material_phasememberAt(grain,ip,el) do source = 1, phase_Nsources(phase) diff --git a/src/constitutive_plastic_disloTungsten.f90 b/src/constitutive_plastic_disloTungsten.f90 index 54c01b912..d9eec28e8 100644 --- a/src/constitutive_plastic_disloTungsten.f90 +++ b/src/constitutive_plastic_disloTungsten.f90 @@ -221,7 +221,7 @@ module function plastic_disloTungsten_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl sizeState = sizeDotState diff --git a/src/constitutive_plastic_dislotwin.f90 b/src/constitutive_plastic_dislotwin.f90 index 62dcdd83e..880c1cb99 100644 --- a/src/constitutive_plastic_dislotwin.f90 +++ b/src/constitutive_plastic_dislotwin.f90 @@ -407,7 +407,7 @@ module function plastic_dislotwin_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl & + size(['f_tw']) * prm%sum_N_tw & + size(['f_tr']) * prm%sum_N_tr diff --git a/src/constitutive_plastic_isotropic.f90 b/src/constitutive_plastic_isotropic.f90 index 2c9028671..1d26fc54c 100644 --- a/src/constitutive_plastic_isotropic.f90 +++ b/src/constitutive_plastic_isotropic.f90 @@ -130,7 +130,7 @@ module function plastic_isotropic_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['xi ','gamma']) sizeState = sizeDotState diff --git a/src/constitutive_plastic_kinehardening.f90 b/src/constitutive_plastic_kinehardening.f90 index 777691242..65c92423c 100644 --- a/src/constitutive_plastic_kinehardening.f90 +++ b/src/constitutive_plastic_kinehardening.f90 @@ -174,7 +174,7 @@ module function plastic_kinehardening_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['crss ','crss_back', 'accshear ']) * prm%sum_N_sl!ToDo: adjust names, ask Philip sizeDeltaState = size(['sense ', 'chi0 ', 'gamma0' ]) * prm%sum_N_sl !ToDo: adjust names sizeState = sizeDotState + sizeDeltaState diff --git a/src/constitutive_plastic_none.f90 b/src/constitutive_plastic_none.f90 index ab5f32d3f..cf2942414 100644 --- a/src/constitutive_plastic_none.f90 +++ b/src/constitutive_plastic_none.f90 @@ -41,7 +41,7 @@ module function plastic_none_init() result(myPlasticity) do p = 1, phases%length phase => phases%get(p) if(.not. myPlasticity(p)) cycle - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs call constitutive_allocateState(plasticState(p),NipcMyPhase,0,0,0) enddo diff --git a/src/constitutive_plastic_nonlocal.f90 b/src/constitutive_plastic_nonlocal.f90 index 86b7df859..a860b21da 100644 --- a/src/constitutive_plastic_nonlocal.f90 +++ b/src/constitutive_plastic_nonlocal.f90 @@ -391,7 +391,7 @@ module function plastic_nonlocal_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs sizeDotState = size([ 'rhoSglEdgePosMobile ','rhoSglEdgeNegMobile ', & 'rhoSglScrewPosMobile ','rhoSglScrewNegMobile ', & 'rhoSglEdgePosImmobile ','rhoSglEdgeNegImmobile ', & @@ -505,7 +505,7 @@ module function plastic_nonlocal_init() result(myPlasticity) enddo allocate(compatibility(2,maxval(param%sum_N_sl),maxval(param%sum_N_sl),nIPneighbors,& - discretization_nIP,discretization_nElem), source=0.0_pReal) + discretization_nIPs,discretization_Nelems), source=0.0_pReal) ! BEGIN DEPRECATED---------------------------------------------------------------------------------- allocate(iRhoU(maxval(param%sum_N_sl),4,Ninstance), source=0) @@ -519,7 +519,7 @@ module function plastic_nonlocal_init() result(myPlasticity) if(.not. myPlasticity(p)) cycle i = i + 1 - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs l = 0 do t = 1,4 do s = 1,param(i)%sum_N_sl @@ -976,7 +976,7 @@ module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, & real(pReal), dimension(3,3), intent(in) :: & Mp !< MandelStress - real(pReal), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems), intent(in) :: & F, & !< elastic deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -1176,7 +1176,7 @@ end subroutine plastic_nonlocal_dotState !--------------------------------------------------------------------------------------------------- function rhoDotFlux(F,Fp,timestep, instance,of,ip,el) - real(pReal), dimension(3,3,homogenization_maxNconstituent,discretization_nIP,discretization_nElem), intent(in) :: & + real(pReal), dimension(3,3,homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems), intent(in) :: & F, & !< elastic deformation gradient Fp !< plastic deformation gradient real(pReal), intent(in) :: & @@ -1416,7 +1416,7 @@ end function rhoDotFlux !-------------------------------------------------------------------------------------------------- module subroutine plastic_nonlocal_updateCompatibility(orientation,instance,i,e) - type(rotation), dimension(1,discretization_nIP,discretization_nElem), intent(in) :: & + type(rotation), dimension(1,discretization_nIPs,discretization_Nelems), intent(in) :: & orientation ! crystal orientation integer, intent(in) :: & instance, & @@ -1632,8 +1632,8 @@ subroutine stateInit(ini,phase,NipcMyPhase,instance) associate(stt => state(instance)) if (ini%random_rho_u > 0.0_pReal) then ! randomly distribute dislocation segments on random slip system and of random type in the volume - do e = 1,discretization_nElem - do i = 1,discretization_nIP + do e = 1,discretization_Nelems + do i = 1,discretization_nIPs if (material_phaseAt(1,e) == phase) volume(material_phasememberAt(1,i,e)) = IPvolume(i,e) enddo enddo diff --git a/src/constitutive_plastic_phenopowerlaw.f90 b/src/constitutive_plastic_phenopowerlaw.f90 index fd1b5fbbb..50c8e835a 100644 --- a/src/constitutive_plastic_phenopowerlaw.f90 +++ b/src/constitutive_plastic_phenopowerlaw.f90 @@ -224,7 +224,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIP + NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['xi_sl ','gamma_sl']) * prm%sum_N_sl & + size(['xi_tw ','gamma_tw']) * prm%sum_N_tw sizeState = sizeDotState diff --git a/src/constitutive_thermal.f90 b/src/constitutive_thermal.f90 index 944e5d72d..a7d5d3259 100644 --- a/src/constitutive_thermal.f90 +++ b/src/constitutive_thermal.f90 @@ -95,7 +95,7 @@ module subroutine constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, homog = material_homogenizationAt(el) instance = thermal_typeInstance(homog) - do grain = 1, homogenization_Nconstituent(homog) + do grain = 1, homogenization_Nconstituents(homog) phase = material_phaseAt(grain,el) constituent = material_phasememberAt(grain,ip,el) do source = 1, phase_Nsources(phase) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index b71bb37f0..c00b076a3 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -135,7 +135,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine crystallite_init - logical, dimension(discretization_nIP,discretization_nElem) :: devNull + logical, dimension(discretization_nIPs,discretization_Nelems) :: devNull integer :: & c, & !< counter in integration point component loop i, & !< counter in integration point loop @@ -162,9 +162,9 @@ subroutine crystallite_init debugCrystallite%ip = config_debug%get_asInt('integrationpoint', defaultVal=1) debugCrystallite%grain = config_debug%get_asInt('grain', defaultVal=1) - cMax = homogenization_maxNconstituent - iMax = discretization_nIP - eMax = discretization_nElem + cMax = homogenization_maxNconstituents + iMax = discretization_nIPs + eMax = discretization_Nelems allocate(crystallite_partitionedF(3,3,cMax,iMax,eMax),source=0.0_pReal) @@ -253,7 +253,7 @@ subroutine crystallite_init ! initialize !$OMP PARALLEL DO PRIVATE(i,c) do e = FEsolving_execElem(1),FEsolving_execElem(2) - do i = FEsolving_execIP(1), FEsolving_execIP(2); do c = 1, homogenization_Nconstituent(material_homogenizationAt(e)) + do i = FEsolving_execIP(1), FEsolving_execIP(2); do c = 1, homogenization_Nconstituents(material_homogenizationAt(e)) crystallite_Fp0(1:3,1:3,c,i,e) = material_orientation0(c,i,e)%asMatrix() ! Fp reflects initial orientation (see 10.1016/j.actamat.2006.01.005) crystallite_Fp0(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) & / math_det33(crystallite_Fp0(1:3,1:3,c,i,e))**(1.0_pReal/3.0_pReal) @@ -279,7 +279,7 @@ subroutine crystallite_init !$OMP PARALLEL DO do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) call constitutive_dependentState(crystallite_partitionedF0(1:3,1:3,c,i,e), & crystallite_partitionedFp0(1:3,1:3,c,i,e), & c,i,e) ! update dependent state variables to be consistent with basic states @@ -307,7 +307,7 @@ end subroutine crystallite_init !-------------------------------------------------------------------------------------------------- function crystallite_stress() - logical, dimension(discretization_nIP,discretization_nElem) :: crystallite_stress + logical, dimension(discretization_nIPs,discretization_Nelems) :: crystallite_stress real(pReal) :: & formerSubStep integer :: & @@ -317,7 +317,7 @@ function crystallite_stress() e, & !< counter in element loop startIP, endIP, & s - logical, dimension(homogenization_maxNconstituent,discretization_nIP,discretization_nElem) :: todo !ToDo: need to set some values to false for different Ngrains + logical, dimension(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems) :: todo !ToDo: need to set some values to false for different Ngrains real(pReal), dimension(:,:,:,:,:), allocatable :: & subLp0,& !< plastic velocity grad at start of crystallite inc subLi0 !< intermediate velocity grad at start of crystallite inc @@ -335,7 +335,7 @@ function crystallite_stress() crystallite_subStep = 0.0_pReal !$OMP PARALLEL DO elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) - do i = FEsolving_execIP(1),FEsolving_execIP(2); do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do i = FEsolving_execIP(1),FEsolving_execIP(2); do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then plasticState (material_phaseAt(c,e))%subState0( :,material_phaseMemberAt(c,i,e)) = & plasticState (material_phaseAt(c,e))%partitionedState0(:,material_phaseMemberAt(c,i,e)) @@ -362,7 +362,7 @@ function crystallite_stress() endIP = startIP else singleRun startIP = 1 - endIP = discretization_nIP + endIP = discretization_nIPs endif singleRun NiterationCrystallite = 0 @@ -376,7 +376,7 @@ function crystallite_stress() !$OMP PARALLEL DO PRIVATE(formerSubStep) elementLooping3: do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) !-------------------------------------------------------------------------------------------------- ! wind forward if (crystallite_converged(c,i,e)) then @@ -472,7 +472,7 @@ subroutine crystallite_initializeRestorationPoints(i,e) c, & !< constituent number s - do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) crystallite_partitionedFp0(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) crystallite_partitionedLp0(1:3,1:3,c,i,e) = crystallite_Lp0(1:3,1:3,c,i,e) crystallite_partitionedFi0(1:3,1:3,c,i,e) = crystallite_Fi0(1:3,1:3,c,i,e) @@ -503,7 +503,7 @@ subroutine crystallite_windForward(i,e) c, & !< constituent number s - do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) crystallite_partitionedF0 (1:3,1:3,c,i,e) = crystallite_partitionedF(1:3,1:3,c,i,e) crystallite_partitionedFp0(1:3,1:3,c,i,e) = crystallite_Fp (1:3,1:3,c,i,e) crystallite_partitionedLp0(1:3,1:3,c,i,e) = crystallite_Lp (1:3,1:3,c,i,e) @@ -536,7 +536,7 @@ subroutine crystallite_restore(i,e,includeL) c, & !< constituent number s - do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) if (includeL) then crystallite_Lp(1:3,1:3,c,i,e) = crystallite_partitionedLp0(1:3,1:3,c,i,e) crystallite_Li(1:3,1:3,c,i,e) = crystallite_partitionedLi0(1:3,1:3,c,i,e) @@ -697,7 +697,7 @@ subroutine crystallite_orientations !$OMP PARALLEL DO do e = FEsolving_execElem(1),FEsolving_execElem(2) do i = FEsolving_execIP(1),FEsolving_execIP(2) - do c = 1,homogenization_Nconstituent(material_homogenizationAt(e)) + do c = 1,homogenization_Nconstituents(material_homogenizationAt(e)) call crystallite_orientation(c,i,e)%fromMatrix(transpose(math_rotationalPart(crystallite_Fe(1:3,1:3,c,i,e)))) enddo; enddo; enddo !$OMP END PARALLEL DO @@ -821,11 +821,11 @@ subroutine crystallite_results real(pReal), allocatable, dimension(:,:,:) :: select_tensors integer :: e,i,c,j - allocate(select_tensors(3,3,count(material_phaseAt==instance)*discretization_nIP)) + allocate(select_tensors(3,3,count(material_phaseAt==instance)*discretization_nIPs)) j=0 do e = 1, size(material_phaseAt,2) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs do c = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains if (material_phaseAt(c,e) == instance) then j = j + 1 @@ -848,11 +848,11 @@ subroutine crystallite_results type(rotation), allocatable, dimension(:) :: select_rotations integer :: e,i,c,j - allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNconstituent*discretization_nIP)) + allocate(select_rotations(count(material_phaseAt==instance)*homogenization_maxNconstituents*discretization_nIPs)) j=0 do e = 1, size(material_phaseAt,2) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs do c = 1, size(material_phaseAt,1) !ToDo: this needs to be changed for varying Ngrains if (material_phaseAt(c,e) == instance) then j = j + 1 diff --git a/src/damage_local.f90 b/src/damage_local.f90 index d36139a99..fefdffce2 100644 --- a/src/damage_local.f90 +++ b/src/damage_local.f90 @@ -143,8 +143,8 @@ subroutine damage_local_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, el call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) - phiDot = phiDot/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) - dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + phiDot = phiDot/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) + dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end subroutine damage_local_getSourceAndItsTangent diff --git a/src/damage_nonlocal.f90 b/src/damage_nonlocal.f90 index f68d778aa..8b494678e 100644 --- a/src/damage_nonlocal.f90 +++ b/src/damage_nonlocal.f90 @@ -110,8 +110,8 @@ subroutine damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, dPhiDot_dPhi = 0.0_pReal call constitutive_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi, ip, el) - phiDot = phiDot/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) - dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + phiDot = phiDot/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) + dPhiDot_dPhi = dPhiDot_dPhi/real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end subroutine damage_nonlocal_getSourceAndItsTangent @@ -132,13 +132,13 @@ function damage_nonlocal_getDiffusion(ip,el) homog = material_homogenizationAt(el) damage_nonlocal_getDiffusion = 0.0_pReal - do grain = 1, homogenization_Nconstituent(homog) + do grain = 1, homogenization_Nconstituents(homog) damage_nonlocal_getDiffusion = damage_nonlocal_getDiffusion + & crystallite_push33ToRef(grain,ip,el,lattice_D(1:3,1:3,material_phaseAt(grain,el))) enddo damage_nonlocal_getDiffusion = & - num%charLength**2*damage_nonlocal_getDiffusion/real(homogenization_Nconstituent(homog),pReal) + num%charLength**2*damage_nonlocal_getDiffusion/real(homogenization_Nconstituents(homog),pReal) end function damage_nonlocal_getDiffusion @@ -156,12 +156,12 @@ real(pReal) function damage_nonlocal_getMobility(ip,el) damage_nonlocal_getMobility = 0.0_pReal - do ipc = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do ipc = 1, homogenization_Nconstituents(material_homogenizationAt(el)) damage_nonlocal_getMobility = damage_nonlocal_getMobility + lattice_M(material_phaseAt(ipc,el)) enddo damage_nonlocal_getMobility = damage_nonlocal_getMobility/& - real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function damage_nonlocal_getMobility diff --git a/src/discretization.f90 b/src/discretization.f90 index e6e53fcf4..0b8925e4a 100644 --- a/src/discretization.f90 +++ b/src/discretization.f90 @@ -11,8 +11,8 @@ module discretization private integer, public, protected :: & - discretization_nIP, & - discretization_nElem + discretization_nIPs, & + discretization_Nelems integer, public, protected, dimension(:), allocatable :: & discretization_materialAt @@ -51,8 +51,8 @@ subroutine discretization_init(materialAt,& print'(/,a)', ' <<<+- discretization init -+>>>'; flush(6) - discretization_nElem = size(materialAt,1) - discretization_nIP = size(IPcoords0,2)/discretization_nElem + discretization_Nelems = size(materialAt,1) + discretization_nIPs = size(IPcoords0,2)/discretization_Nelems discretization_materialAt = materialAt diff --git a/src/homogenization.f90 b/src/homogenization.f90 index c717c24fd..02f1d3bf0 100644 --- a/src/homogenization.f90 +++ b/src/homogenization.f90 @@ -158,7 +158,7 @@ subroutine homogenization_init debugHomog%grain = config_debug%get_asInt('grain',defaultVal = 1) if (debugHomog%grain < 1 & - .or. debugHomog%grain > homogenization_Nconstituent(material_homogenizationAt(debugHomog%element))) & + .or. debugHomog%grain > homogenization_Nconstituents(material_homogenizationAt(debugHomog%element))) & call IO_error(602,ext_msg='constituent', el=debugHomog%element, g=debugHomog%grain) @@ -181,10 +181,10 @@ subroutine homogenization_init !-------------------------------------------------------------------------------------------------- ! allocate and initialize global variables - allocate(homogenization_dPdF(3,3,3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) - homogenization_F0 = spread(spread(math_I3,3,discretization_nIP),4,discretization_nElem) ! initialize to identity + allocate(homogenization_dPdF(3,3,3,3,discretization_nIPs,discretization_Nelems), source=0.0_pReal) + homogenization_F0 = spread(spread(math_I3,3,discretization_nIPs),4,discretization_Nelems) ! initialize to identity homogenization_F = homogenization_F0 ! initialize to identity - allocate(homogenization_P(3,3,discretization_nIP,discretization_nElem), source=0.0_pReal) + allocate(homogenization_P(3,3,discretization_nIPs,discretization_Nelems), source=0.0_pReal) print'(/,a)', ' <<<+- homogenization init -+>>>'; flush(IO_STDOUT) @@ -213,13 +213,13 @@ subroutine materialpoint_stressAndItsTangent(dt) i, & !< integration point number e, & !< element number myNgrains - real(pReal), dimension(discretization_nIP,discretization_nElem) :: & + real(pReal), dimension(discretization_nIPs,discretization_Nelems) :: & subFrac, & subStep - logical, dimension(discretization_nIP,discretization_nElem) :: & + logical, dimension(discretization_nIPs,discretization_Nelems) :: & requested, & converged - logical, dimension(2,discretization_nIP,discretization_nElem) :: & + logical, dimension(2,discretization_nIPs,discretization_Nelems) :: & doneAndHappy @@ -257,7 +257,7 @@ subroutine materialpoint_stressAndItsTangent(dt) !$OMP PARALLEL DO elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) - myNgrains = homogenization_Nconstituent(material_homogenizationAt(e)) + myNgrains = homogenization_Nconstituents(material_homogenizationAt(e)) IpLooping1: do i = FEsolving_execIP(1),FEsolving_execIP(2) if (converged(i,e)) then @@ -327,7 +327,7 @@ subroutine materialpoint_stressAndItsTangent(dt) ! deformation partitioning !$OMP PARALLEL DO PRIVATE(myNgrains) elementLooping2: do e = FEsolving_execElem(1),FEsolving_execElem(2) - myNgrains = homogenization_Nconstituent(material_homogenizationAt(e)) + myNgrains = homogenization_Nconstituents(material_homogenizationAt(e)) IpLooping2: do i = FEsolving_execIP(1),FEsolving_execIP(2) if(requested(i,e) .and. .not. doneAndHappy(1,i,e)) then ! requested but not yet done call partitionDeformation(homogenization_F0(1:3,1:3,i,e) & @@ -408,12 +408,12 @@ subroutine partitionDeformation(subF,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization call mech_isostrain_partitionDeformation(& - crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & subF) case (HOMOGENIZATION_RGC_ID) chosenHomogenization call mech_RGC_partitionDeformation(& - crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & subF,& ip, & el) @@ -437,19 +437,19 @@ function updateState(subdt,subF,ip,el) el !< element number integer :: c logical, dimension(2) :: updateState - real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituent(material_homogenizationAt(el))) + real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituents(material_homogenizationAt(el))) updateState = .true. chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - do c=1,homogenization_Nconstituent(material_homogenizationAt(el)) + do c=1,homogenization_Nconstituents(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo updateState = & updateState .and. & - mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & - crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & - crystallite_partitionedF0(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el),& + mech_RGC_updateState(crystallite_P(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & + crystallite_partitionedF0(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el),& subF,& subdt, & dPdFs, & @@ -487,7 +487,7 @@ subroutine averageStressAndItsTangent(ip,el) ip, & !< integration point el !< element number integer :: c - real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituent(material_homogenizationAt(el))) + real(pReal) :: dPdFs(3,3,3,3,homogenization_Nconstituents(material_homogenizationAt(el))) chosenHomogenization: select case(homogenization_type(material_homogenizationAt(el))) @@ -496,24 +496,24 @@ subroutine averageStressAndItsTangent(ip,el) homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el) = crystallite_stressTangent(1,ip,el) case (HOMOGENIZATION_ISOSTRAIN_ID) chosenHomogenization - do c = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do c = 1, homogenization_Nconstituents(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_isostrain_averageStressAndItsTangent(& homogenization_P(1:3,1:3,ip,el), & homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el),& - crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & + crystallite_P(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) case (HOMOGENIZATION_RGC_ID) chosenHomogenization - do c = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do c = 1, homogenization_Nconstituents(material_homogenizationAt(el)) dPdFs(:,:,:,:,c) = crystallite_stressTangent(c,ip,el) enddo call mech_RGC_averageStressAndItsTangent(& homogenization_P(1:3,1:3,ip,el), & homogenization_dPdF(1:3,1:3,1:3,1:3,ip,el),& - crystallite_P(1:3,1:3,1:homogenization_Nconstituent(material_homogenizationAt(el)),ip,el), & + crystallite_P(1:3,1:3,1:homogenization_Nconstituents(material_homogenizationAt(el)),ip,el), & dPdFs, & homogenization_typeInstance(material_homogenizationAt(el))) end select chosenHomogenization @@ -539,10 +539,10 @@ subroutine homogenization_results group = trim(group_base)//'/generic' call results_closeGroup(results_addGroup(group)) - !temp = reshape(homogenization_F,[3,3,discretization_nIP*discretization_nElem]) + !temp = reshape(homogenization_F,[3,3,discretization_nIPs*discretization_Nelems]) !call results_writeDataset(group,temp,'F',& ! 'deformation gradient','1') - !temp = reshape(homogenization_P,[3,3,discretization_nIP*discretization_nElem]) + !temp = reshape(homogenization_P,[3,3,discretization_nIPs*discretization_Nelems]) !call results_writeDataset(group,temp,'P',& ! '1st Piola-Kirchhoff stress','Pa') diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index 6bbbbbaa9..1a0f7a5d3 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -164,7 +164,7 @@ module subroutine mech_RGC_init(num_homogMech) #endif prm%N_constituents = homogMech%get_asInts('cluster_size',requiredSize=3) - if (homogenization_Nconstituent(h) /= product(prm%N_constituents)) & + if (homogenization_Nconstituents(h) /= product(prm%N_constituents)) & call IO_error(211,ext_msg='N_constituents (mech_RGC)') prm%xi_alpha = homogMech%get_asFloat('xi_alpha') diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index 5a03e1204..7f3724ae1 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -51,7 +51,7 @@ module subroutine mech_isostrain_init homogMech => homog%get('mech') associate(prm => param(homogenization_typeInstance(h))) - prm%N_constituents = homogenization_Nconstituent(h) + prm%N_constituents = homogenization_Nconstituents(h) select case(homogMech%get_asString('mapping',defaultVal = 'sum')) case ('sum') prm%mapping = parallel_ID diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index d64154dcc..80949010e 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -26,7 +26,7 @@ module subroutine mech_none_init do h = 1, size(homogenization_type) if(homogenization_type(h) /= HOMOGENIZATION_NONE_ID) cycle - if(homogenization_Nconstituent(h) /= 1) & + if(homogenization_Nconstituents(h) /= 1) & call IO_error(211,ext_msg='N_constituents (mech_none)') NofMyHomog = count(material_homogenizationAt == h) diff --git a/src/material.f90 b/src/material.f90 index 58a5c4080..a0c81a519 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -64,10 +64,10 @@ module material homogenization_type !< type of each homogenization integer, public, protected :: & - homogenization_maxNconstituent !< max number of grains in any USED homogenization + homogenization_maxNconstituents !< max number of grains in any USED homogenization integer, dimension(:), allocatable, public, protected :: & - homogenization_Nconstituent, & !< number of grains in each homogenization + homogenization_Nconstituents, & !< number of grains in each 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 @@ -183,7 +183,7 @@ subroutine material_init(restart) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! BEGIN DEPRECATED - allocate(mappingHomogenizationConst( discretization_nIP,discretization_nElem),source=1) + allocate(mappingHomogenizationConst( discretization_nIPs,discretization_Nelems),source=1) ! hack needed to initialize field values used during constitutive initialization do myHomog = 1, size(material_name_homogenization) @@ -312,29 +312,29 @@ subroutine material_parseMaterial material_name_phase = getKeys(phases) material_name_homogenization = getKeys(homogenizations) - allocate(homogenization_Nconstituent(homogenizations%length)) + allocate(homogenization_Nconstituents(homogenizations%length)) do h=1, homogenizations%length homogenization => homogenizations%get(h) - homogenization_Nconstituent(h) = homogenization%get_asInt('N_constituents') + homogenization_Nconstituents(h) = homogenization%get_asInt('N_constituents') enddo - homogenization_maxNconstituent = maxval(homogenization_Nconstituent) + homogenization_maxNconstituents = maxval(homogenization_Nconstituents) allocate(counterPhase(phases%length),source=0) allocate(counterHomogenization(homogenizations%length),source=0) - allocate(material_homogenizationAt(discretization_nElem),source=0) - allocate(material_homogenizationMemberAt(discretization_nIP,discretization_nElem),source=0) - allocate(material_phaseAt(homogenization_maxNconstituent,discretization_nElem),source=0) - allocate(material_phaseMemberAt(homogenization_maxNconstituent,discretization_nIP,discretization_nElem),source=0) + allocate(material_homogenizationAt(discretization_Nelems),source=0) + allocate(material_homogenizationMemberAt(discretization_nIPs,discretization_Nelems),source=0) + allocate(material_phaseAt(homogenization_maxNconstituents,discretization_Nelems),source=0) + allocate(material_phaseMemberAt(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems),source=0) - allocate(material_orientation0(homogenization_maxNconstituent,discretization_nIP,discretization_nElem)) + allocate(material_orientation0(homogenization_maxNconstituents,discretization_nIPs,discretization_Nelems)) - do e = 1, discretization_nElem + do e = 1, discretization_Nelems material => materials%get(discretization_materialAt(e)) constituents => material%get('constituents') material_homogenizationAt(e) = homogenizations%getIndex(material%get_asString('homogenization')) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs counterHomogenization(material_homogenizationAt(e)) = counterHomogenization(material_homogenizationAt(e)) + 1 material_homogenizationMemberAt(i,e) = counterHomogenization(material_homogenizationAt(e)) enddo @@ -345,7 +345,7 @@ subroutine material_parseMaterial frac = frac + constituent%get_asFloat('fraction') material_phaseAt(c,e) = phases%getIndex(constituent%get_asString('phase')) - do i = 1, discretization_nIP + do i = 1, discretization_nIPs counterPhase(material_phaseAt(c,e)) = counterPhase(material_phaseAt(c,e)) + 1 material_phaseMemberAt(c,i,e) = counterPhase(material_phaseAt(c,e)) diff --git a/src/source_damage_anisoBrittle.f90 b/src/source_damage_anisoBrittle.f90 index 1954e6ea3..faeb8bd87 100644 --- a/src/source_damage_anisoBrittle.f90 +++ b/src/source_damage_anisoBrittle.f90 @@ -100,7 +100,7 @@ module function source_damage_anisoBrittle_init(source_length) result(mySources) if (any(prm%g_crit < 0.0_pReal)) extmsg = trim(extmsg)//' g_crit' if (any(prm%s_crit < 0.0_pReal)) extmsg = trim(extmsg)//' s_crit' - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol' diff --git a/src/source_damage_anisoDuctile.f90 b/src/source_damage_anisoDuctile.f90 index 3b8eb2428..fd3fa38ed 100644 --- a/src/source_damage_anisoDuctile.f90 +++ b/src/source_damage_anisoDuctile.f90 @@ -84,7 +84,7 @@ module function source_damage_anisoDuctile_init(source_length) result(mySources) if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' if (any(prm%gamma_crit < 0.0_pReal)) extmsg = trim(extmsg)//' gamma_crit' - NipcMyPhase=count(material_phaseAt==p) * discretization_nIP + NipcMyPhase=count(material_phaseAt==p) * discretization_nIPs call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('anisoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_atol' diff --git a/src/source_damage_isoBrittle.f90 b/src/source_damage_isoBrittle.f90 index 714e71ef1..c56407c3d 100644 --- a/src/source_damage_isoBrittle.f90 +++ b/src/source_damage_isoBrittle.f90 @@ -73,7 +73,7 @@ module function source_damage_isoBrittle_init(source_length) result(mySources) ! sanity checks if (prm%W_crit <= 0.0_pReal) extmsg = trim(extmsg)//' W_crit' - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,1) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('isoBrittle_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isobrittle_atol' diff --git a/src/source_damage_isoDuctile.f90 b/src/source_damage_isoDuctile.f90 index d958aed6a..0d068b6e4 100644 --- a/src/source_damage_isoDuctile.f90 +++ b/src/source_damage_isoDuctile.f90 @@ -77,7 +77,7 @@ module function source_damage_isoDuctile_init(source_length) result(mySources) if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' if (prm%gamma_crit <= 0.0_pReal) extmsg = trim(extmsg)//' gamma_crit' - NipcMyPhase=count(material_phaseAt==p) * discretization_nIP + NipcMyPhase=count(material_phaseAt==p) * discretization_nIPs call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('isoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isoductile_atol' diff --git a/src/source_thermal_dissipation.f90 b/src/source_thermal_dissipation.f90 index 5cc740424..734451a72 100644 --- a/src/source_thermal_dissipation.f90 +++ b/src/source_thermal_dissipation.f90 @@ -61,7 +61,7 @@ module function source_thermal_dissipation_init(source_length) result(mySources) src => sources%get(sourceOffset) prm%kappa = src%get_asFloat('kappa') - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,0,0,0) end associate diff --git a/src/source_thermal_externalheat.f90 b/src/source_thermal_externalheat.f90 index 2eeeb47df..cbc1fa69d 100644 --- a/src/source_thermal_externalheat.f90 +++ b/src/source_thermal_externalheat.f90 @@ -69,7 +69,7 @@ module function source_thermal_externalheat_init(source_length) result(mySources prm%f_T = src%get_asFloats('f_T',requiredSize = size(prm%t_n)) - NipcMyPhase = count(material_phaseAt==p) * discretization_nIP + NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) end associate diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index 07dbdebff..189a7131a 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -145,8 +145,8 @@ subroutine thermal_adiabatic_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) homog = material_homogenizationAt(el) call constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, crystallite_S, crystallite_Lp, ip, el) - Tdot = Tdot/real(homogenization_Nconstituent(homog),pReal) - dTdot_dT = dTdot_dT/real(homogenization_Nconstituent(homog),pReal) + Tdot = Tdot/real(homogenization_Nconstituents(homog),pReal) + dTdot_dT = dTdot_dT/real(homogenization_Nconstituents(homog),pReal) end subroutine thermal_adiabatic_getSourceAndItsTangent @@ -167,13 +167,13 @@ function thermal_adiabatic_getSpecificHeat(ip,el) thermal_adiabatic_getSpecificHeat = 0.0_pReal - do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_adiabatic_getSpecificHeat = thermal_adiabatic_getSpecificHeat & + lattice_c_p(material_phaseAt(grain,el)) enddo thermal_adiabatic_getSpecificHeat = thermal_adiabatic_getSpecificHeat & - / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_adiabatic_getSpecificHeat @@ -193,13 +193,13 @@ function thermal_adiabatic_getMassDensity(ip,el) thermal_adiabatic_getMassDensity = 0.0_pReal - do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_adiabatic_getMassDensity = thermal_adiabatic_getMassDensity & + lattice_rho(material_phaseAt(grain,el)) enddo thermal_adiabatic_getMassDensity = thermal_adiabatic_getMassDensity & - / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_adiabatic_getMassDensity diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index 616965df0..ead06306a 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -104,8 +104,8 @@ subroutine thermal_conduction_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el) homog = material_homogenizationAt(el) call constitutive_thermal_getRateAndItsTangents(TDot, dTDot_dT, T, crystallite_S,crystallite_Lp ,ip, el) - Tdot = Tdot/real(homogenization_Nconstituent(homog),pReal) - dTdot_dT = dTdot_dT/real(homogenization_Nconstituent(homog),pReal) + Tdot = Tdot/real(homogenization_Nconstituents(homog),pReal) + dTdot_dT = dTdot_dT/real(homogenization_Nconstituents(homog),pReal) end subroutine thermal_conduction_getSourceAndItsTangent @@ -125,13 +125,13 @@ function thermal_conduction_getConductivity(ip,el) thermal_conduction_getConductivity = 0.0_pReal - do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_conduction_getConductivity = thermal_conduction_getConductivity + & crystallite_push33ToRef(grain,ip,el,lattice_K(:,:,material_phaseAt(grain,el))) enddo thermal_conduction_getConductivity = thermal_conduction_getConductivity & - / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_conduction_getConductivity @@ -151,13 +151,13 @@ function thermal_conduction_getSpecificHeat(ip,el) thermal_conduction_getSpecificHeat = 0.0_pReal - do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_conduction_getSpecificHeat = thermal_conduction_getSpecificHeat & + lattice_c_p(material_phaseAt(grain,el)) enddo thermal_conduction_getSpecificHeat = thermal_conduction_getSpecificHeat & - / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_conduction_getSpecificHeat @@ -178,13 +178,13 @@ function thermal_conduction_getMassDensity(ip,el) thermal_conduction_getMassDensity = 0.0_pReal - do grain = 1, homogenization_Nconstituent(material_homogenizationAt(el)) + do grain = 1, homogenization_Nconstituents(material_homogenizationAt(el)) thermal_conduction_getMassDensity = thermal_conduction_getMassDensity & + lattice_rho(material_phaseAt(grain,el)) enddo thermal_conduction_getMassDensity = thermal_conduction_getMassDensity & - / real(homogenization_Nconstituent(material_homogenizationAt(el)),pReal) + / real(homogenization_Nconstituents(material_homogenizationAt(el)),pReal) end function thermal_conduction_getMassDensity From 568d1a020cd64dd4ab605d4e3849d666b53885a5 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Tue, 27 Oct 2020 21:33:30 +0100 Subject: [PATCH 19/36] consistent names for counting variables --- src/constitutive.f90 | 16 ++--- src/constitutive_plastic_disloTungsten.f90 | 30 ++++----- src/constitutive_plastic_dislotwin.f90 | 46 +++++++------- src/constitutive_plastic_isotropic.f90 | 20 +++--- src/constitutive_plastic_kinehardening.f90 | 24 ++++---- src/constitutive_plastic_none.f90 | 14 ++--- src/constitutive_plastic_nonlocal.f90 | 72 +++++++++++----------- src/constitutive_plastic_phenopowerlaw.f90 | 24 ++++---- src/damage_local.f90 | 14 ++--- src/damage_none.f90 | 10 +-- src/damage_nonlocal.f90 | 14 ++--- src/homogenization_mech_RGC.f90 | 36 +++++------ src/homogenization_mech_isostrain.f90 | 20 +++--- src/homogenization_mech_none.f90 | 16 ++--- src/kinematics_cleavage_opening.f90 | 12 ++-- src/kinematics_slipplane_opening.f90 | 12 ++-- src/kinematics_thermal_expansion.f90 | 10 +-- src/source_damage_anisoBrittle.f90 | 16 ++--- src/source_damage_anisoDuctile.f90 | 16 ++--- src/source_damage_isoBrittle.f90 | 16 ++--- src/source_damage_isoDuctile.f90 | 16 ++--- src/source_thermal_dissipation.f90 | 18 +++--- src/source_thermal_externalheat.f90 | 16 ++--- src/thermal_adiabatic.f90 | 18 +++--- src/thermal_conduction.f90 | 18 +++--- src/thermal_isothermal.f90 | 10 +-- 26 files changed, 267 insertions(+), 267 deletions(-) diff --git a/src/constitutive.f90 b/src/constitutive.f90 index c7658b77f..358937e4b 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -905,12 +905,12 @@ end function constitutive_deltaState !> @brief Allocate the components of the state structure for a given phase !-------------------------------------------------------------------------------------------------- subroutine constitutive_allocateState(state, & - NipcMyPhase,sizeState,sizeDotState,sizeDeltaState) + Nconstituents,sizeState,sizeDotState,sizeDeltaState) class(tState), intent(out) :: & state integer, intent(in) :: & - NipcMyPhase, & + Nconstituents, & sizeState, & sizeDotState, & sizeDeltaState @@ -921,14 +921,14 @@ subroutine constitutive_allocateState(state, & state%offsetDeltaState = sizeState-sizeDeltaState ! deltaState occupies latter part of state by definition allocate(state%atol (sizeState), source=0.0_pReal) - allocate(state%state0 (sizeState,NipcMyPhase), source=0.0_pReal) - allocate(state%partitionedState0(sizeState,NipcMyPhase), source=0.0_pReal) - allocate(state%subState0 (sizeState,NipcMyPhase), source=0.0_pReal) - allocate(state%state (sizeState,NipcMyPhase), source=0.0_pReal) + allocate(state%state0 (sizeState,Nconstituents), source=0.0_pReal) + allocate(state%partitionedState0(sizeState,Nconstituents), source=0.0_pReal) + allocate(state%subState0 (sizeState,Nconstituents), source=0.0_pReal) + allocate(state%state (sizeState,Nconstituents), source=0.0_pReal) - allocate(state%dotState (sizeDotState,NipcMyPhase), source=0.0_pReal) + allocate(state%dotState (sizeDotState,Nconstituents), source=0.0_pReal) - allocate(state%deltaState(sizeDeltaState,NipcMyPhase), source=0.0_pReal) + allocate(state%deltaState(sizeDeltaState,Nconstituents), source=0.0_pReal) end subroutine constitutive_allocateState diff --git a/src/constitutive_plastic_disloTungsten.f90 b/src/constitutive_plastic_disloTungsten.f90 index d9eec28e8..aef3bef89 100644 --- a/src/constitutive_plastic_disloTungsten.f90 +++ b/src/constitutive_plastic_disloTungsten.f90 @@ -78,9 +78,9 @@ module function plastic_disloTungsten_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -99,17 +99,17 @@ module function plastic_disloTungsten_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_dislotungsten init -+>>>' myPlasticity = plastic_active('disloTungsten') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return print*, 'Cereceda et al., International Journal of Plasticity 78:242–256, 2016' print*, 'https://dx.doi.org/10.1016/j.ijplas.2015.09.002' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) - allocate(dependentState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) + allocate(dependentState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -221,18 +221,18 @@ module function plastic_disloTungsten_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization startIndex = 1 endIndex = prm%sum_N_sl stt%rho_mob => plasticState(p)%state(startIndex:endIndex,:) - stt%rho_mob = spread(rho_mob_0,2,NipcMyPhase) + stt%rho_mob = spread(rho_mob_0,2,Nconstituents) dot%rho_mob => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) if (any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_rho' @@ -240,7 +240,7 @@ module function plastic_disloTungsten_init() result(myPlasticity) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl stt%rho_dip => plasticState(p)%state(startIndex:endIndex,:) - stt%rho_dip = spread(rho_dip_0,2,NipcMyPhase) + stt%rho_dip = spread(rho_dip_0,2,Nconstituents) dot%rho_dip => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) @@ -252,8 +252,8 @@ module function plastic_disloTungsten_init() result(myPlasticity) ! global alias plasticState(p)%slipRate => plasticState(p)%dotState(startIndex:endIndex,:) - allocate(dst%Lambda_sl(prm%sum_N_sl,NipcMyPhase), source=0.0_pReal) - allocate(dst%threshold_stress(prm%sum_N_sl,NipcMyPhase), source=0.0_pReal) + allocate(dst%Lambda_sl(prm%sum_N_sl,Nconstituents), source=0.0_pReal) + allocate(dst%threshold_stress(prm%sum_N_sl,Nconstituents), source=0.0_pReal) plasticState(p)%state0 = plasticState(p)%state ! ToDo: this could be done centrally diff --git a/src/constitutive_plastic_dislotwin.f90 b/src/constitutive_plastic_dislotwin.f90 index 880c1cb99..30527e301 100644 --- a/src/constitutive_plastic_dislotwin.f90 +++ b/src/constitutive_plastic_dislotwin.f90 @@ -126,9 +126,9 @@ module function plastic_dislotwin_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -146,9 +146,9 @@ module function plastic_dislotwin_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_dislotwin init -+>>>' myPlasticity = plastic_active('dislotwin') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return print*, 'Ma and Roters, Acta Materialia 52(12):3603–3612, 2004' print*, 'https://doi.org/10.1016/j.actamat.2004.04.012'//IO_EOL @@ -159,10 +159,10 @@ module function plastic_dislotwin_init() result(myPlasticity) print*, 'Wong et al., Acta Materialia 118:140–151, 2016' print*, 'https://doi.org/10.1016/j.actamat.2016.07.032' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) - allocate(dependentState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) + allocate(dependentState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -407,21 +407,21 @@ module function plastic_dislotwin_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl & + size(['f_tw']) * prm%sum_N_tw & + size(['f_tr']) * prm%sum_N_tr sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! locally defined state aliases and initialization of state0 and atol startIndex = 1 endIndex = prm%sum_N_sl stt%rho_mob=>plasticState(p)%state(startIndex:endIndex,:) - stt%rho_mob= spread(rho_mob_0,2,NipcMyPhase) + stt%rho_mob= spread(rho_mob_0,2,Nconstituents) dot%rho_mob=>plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) if (any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_rho' @@ -429,7 +429,7 @@ module function plastic_dislotwin_init() result(myPlasticity) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_sl stt%rho_dip=>plasticState(p)%state(startIndex:endIndex,:) - stt%rho_dip= spread(rho_dip_0,2,NipcMyPhase) + stt%rho_dip= spread(rho_dip_0,2,Nconstituents) dot%rho_dip=>plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_rho',defaultVal=1.0_pReal) @@ -455,18 +455,18 @@ module function plastic_dislotwin_init() result(myPlasticity) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('f_trans',defaultVal=1.0e-6_pReal) if (any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' f_trans' - allocate(dst%Lambda_sl (prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_pass (prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) + allocate(dst%Lambda_sl (prm%sum_N_sl,Nconstituents),source=0.0_pReal) + allocate(dst%tau_pass (prm%sum_N_sl,Nconstituents),source=0.0_pReal) - allocate(dst%Lambda_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_hat_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_r_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) - allocate(dst%V_tw (prm%sum_N_tw,NipcMyPhase),source=0.0_pReal) + allocate(dst%Lambda_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) + allocate(dst%tau_hat_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) + allocate(dst%tau_r_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) + allocate(dst%V_tw (prm%sum_N_tw,Nconstituents),source=0.0_pReal) - allocate(dst%Lambda_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_hat_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_r_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) - allocate(dst%V_tr (prm%sum_N_tr,NipcMyPhase),source=0.0_pReal) + allocate(dst%Lambda_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) + allocate(dst%tau_hat_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) + allocate(dst%tau_r_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) + allocate(dst%V_tr (prm%sum_N_tr,Nconstituents),source=0.0_pReal) plasticState(p)%state0 = plasticState(p)%state ! ToDo: this could be done centrally diff --git a/src/constitutive_plastic_isotropic.f90 b/src/constitutive_plastic_isotropic.f90 index 1d26fc54c..db418663f 100644 --- a/src/constitutive_plastic_isotropic.f90 +++ b/src/constitutive_plastic_isotropic.f90 @@ -53,10 +53,10 @@ module function plastic_isotropic_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, & i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState real(pReal) :: & xi_0 !< initial critical stress @@ -70,16 +70,16 @@ module function plastic_isotropic_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_isotropic init -+>>>' myPlasticity = plastic_active('isotropic') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return print*, 'Maiti and Eisenlohr, Scripta Materialia 145:37–40, 2018' print*, 'https://doi.org/10.1016/j.scriptamat.2017.09.047' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -130,11 +130,11 @@ module function plastic_isotropic_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['xi ','gamma']) sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization diff --git a/src/constitutive_plastic_kinehardening.f90 b/src/constitutive_plastic_kinehardening.f90 index 65c92423c..3faf7dc41 100644 --- a/src/constitutive_plastic_kinehardening.f90 +++ b/src/constitutive_plastic_kinehardening.f90 @@ -62,9 +62,9 @@ module function plastic_kinehardening_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, o, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDeltaState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -82,14 +82,14 @@ module function plastic_kinehardening_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_kinehardening init -+>>>' myPlasticity = plastic_active('kinehardening') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) - allocate(deltaState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) + allocate(deltaState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -174,19 +174,19 @@ module function plastic_kinehardening_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['crss ','crss_back', 'accshear ']) * prm%sum_N_sl!ToDo: adjust names, ask Philip sizeDeltaState = size(['sense ', 'chi0 ', 'gamma0' ]) * prm%sum_N_sl !ToDo: adjust names sizeState = sizeDotState + sizeDeltaState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,sizeDeltaState) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization startIndex = 1 endIndex = prm%sum_N_sl stt%crss => plasticState(p)%state (startIndex:endIndex,:) - stt%crss = spread(xi_0, 2, NipcMyPhase) + stt%crss = spread(xi_0, 2, Nconstituents) dot%crss => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) if(any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' diff --git a/src/constitutive_plastic_none.f90 b/src/constitutive_plastic_none.f90 index cf2942414..923ae354a 100644 --- a/src/constitutive_plastic_none.f90 +++ b/src/constitutive_plastic_none.f90 @@ -16,9 +16,9 @@ module function plastic_none_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, & - NipcMyPhase + Nconstituents class(tNode), pointer :: & phases, & phase, & @@ -34,15 +34,15 @@ module function plastic_none_init() result(myPlasticity) if(pl%get_asString('type') == 'none') myPlasticity(p) = .true. enddo - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return do p = 1, phases%length phase => phases%get(p) if(.not. myPlasticity(p)) cycle - NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs - call constitutive_allocateState(plasticState(p),NipcMyPhase,0,0,0) + Nconstituents = count(material_phaseAt == p) * discretization_nIPs + call constitutive_allocateState(plasticState(p),Nconstituents,0,0,0) enddo end function plastic_none_init diff --git a/src/constitutive_plastic_nonlocal.f90 b/src/constitutive_plastic_nonlocal.f90 index a860b21da..472415ba9 100644 --- a/src/constitutive_plastic_nonlocal.f90 +++ b/src/constitutive_plastic_nonlocal.f90 @@ -153,7 +153,7 @@ submodule(constitutive:constitutive_plastic) plastic_nonlocal state, & state0 - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) type(tNonlocalMicrostructure), dimension(:), allocatable :: microstructure @@ -168,9 +168,9 @@ module function plastic_nonlocal_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, sizeDependentState, sizeDeltaState, & s1, s2, & s, t, l @@ -188,9 +188,9 @@ module function plastic_nonlocal_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_nonlocal init -+>>>' myPlasticity = plastic_active('nonlocal') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) then + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) then call geometry_plastic_nonlocal_disable return endif @@ -201,12 +201,12 @@ module function plastic_nonlocal_init() result(myPlasticity) print*, 'Kords, Dissertation RWTH Aachen, 2014' print*, 'http://publications.rwth-aachen.de/record/229993' - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(state0(Ninstance)) - allocate(dotState(Ninstance)) - allocate(deltaState(Ninstance)) - allocate(microstructure(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(state0(Ninstances)) + allocate(dotState(Ninstances)) + allocate(deltaState(Ninstances)) + allocate(microstructure(Ninstances)) phases => config_material%get('phase') i = 0 @@ -391,7 +391,7 @@ module function plastic_nonlocal_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs + Nconstituents = count(material_phaseAt==p) * discretization_nIPs sizeDotState = size([ 'rhoSglEdgePosMobile ','rhoSglEdgeNegMobile ', & 'rhoSglScrewPosMobile ','rhoSglScrewNegMobile ', & 'rhoSglEdgePosImmobile ','rhoSglEdgeNegImmobile ', & @@ -405,7 +405,7 @@ module function plastic_nonlocal_init() result(myPlasticity) 'maxDipoleHeightEdge ','maxDipoleHeightScrew' ]) * prm%sum_N_sl !< other dependent state variables that are not updated by microstructure sizeDeltaState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,sizeDeltaState) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,sizeDeltaState) plasticState(p)%nonlocal = pl%get_asBool('nonlocal') if(plasticState(p)%nonlocal .and. .not. allocated(IPneighborhood)) & @@ -476,26 +476,26 @@ module function plastic_nonlocal_init() result(myPlasticity) dot%rho_dip_scr => plasticState(p)%dotState (9*prm%sum_N_sl+1:10*prm%sum_N_sl,:) del%rho_dip_scr => plasticState(p)%deltaState (9*prm%sum_N_sl+1:10*prm%sum_N_sl,:) - stt%gamma => plasticState(p)%state (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) - dot%gamma => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) - del%gamma => plasticState(p)%deltaState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) + stt%gamma => plasticState(p)%state (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) + dot%gamma => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) + del%gamma => plasticState(p)%deltaState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) plasticState(p)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl ) = pl%get_asFloat('atol_gamma', defaultVal = 1.0e-2_pReal) if(any(plasticState(p)%atol(10*prm%sum_N_sl+1:11*prm%sum_N_sl) < 0.0_pReal)) & extmsg = trim(extmsg)//' atol_gamma' - plasticState(p)%slipRate => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:NipcMyPhase) + plasticState(p)%slipRate => plasticState(p)%dotState (10*prm%sum_N_sl + 1:11*prm%sum_N_sl,1:Nconstituents) - stt%rho_forest => plasticState(p)%state (11*prm%sum_N_sl + 1:12*prm%sum_N_sl,1:NipcMyPhase) - stt%v => plasticState(p)%state (12*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:NipcMyPhase) - stt%v_edg_pos => plasticState(p)%state (12*prm%sum_N_sl + 1:13*prm%sum_N_sl,1:NipcMyPhase) - stt%v_edg_neg => plasticState(p)%state (13*prm%sum_N_sl + 1:14*prm%sum_N_sl,1:NipcMyPhase) - stt%v_scr_pos => plasticState(p)%state (14*prm%sum_N_sl + 1:15*prm%sum_N_sl,1:NipcMyPhase) - stt%v_scr_neg => plasticState(p)%state (15*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:NipcMyPhase) + stt%rho_forest => plasticState(p)%state (11*prm%sum_N_sl + 1:12*prm%sum_N_sl,1:Nconstituents) + stt%v => plasticState(p)%state (12*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:Nconstituents) + stt%v_edg_pos => plasticState(p)%state (12*prm%sum_N_sl + 1:13*prm%sum_N_sl,1:Nconstituents) + stt%v_edg_neg => plasticState(p)%state (13*prm%sum_N_sl + 1:14*prm%sum_N_sl,1:Nconstituents) + stt%v_scr_pos => plasticState(p)%state (14*prm%sum_N_sl + 1:15*prm%sum_N_sl,1:Nconstituents) + stt%v_scr_neg => plasticState(p)%state (15*prm%sum_N_sl + 1:16*prm%sum_N_sl,1:Nconstituents) - allocate(dst%tau_pass(prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) - allocate(dst%tau_back(prm%sum_N_sl,NipcMyPhase),source=0.0_pReal) + allocate(dst%tau_pass(prm%sum_N_sl,Nconstituents),source=0.0_pReal) + allocate(dst%tau_back(prm%sum_N_sl,Nconstituents),source=0.0_pReal) end associate - if (NipcMyPhase > 0) call stateInit(ini,p,NipcMyPhase,i) + if (Nconstituents > 0) call stateInit(ini,p,Nconstituents,i) plasticState(p)%state0 = plasticState(p)%state !-------------------------------------------------------------------------------------------------- @@ -508,9 +508,9 @@ module function plastic_nonlocal_init() result(myPlasticity) discretization_nIPs,discretization_Nelems), source=0.0_pReal) ! BEGIN DEPRECATED---------------------------------------------------------------------------------- - allocate(iRhoU(maxval(param%sum_N_sl),4,Ninstance), source=0) - allocate(iV(maxval(param%sum_N_sl),4,Ninstance), source=0) - allocate(iD(maxval(param%sum_N_sl),2,Ninstance), source=0) + allocate(iRhoU(maxval(param%sum_N_sl),4,Ninstances), source=0) + allocate(iV(maxval(param%sum_N_sl),4,Ninstances), source=0) + allocate(iD(maxval(param%sum_N_sl),2,Ninstances), source=0) i = 0 do p = 1, phases%length @@ -519,7 +519,7 @@ module function plastic_nonlocal_init() result(myPlasticity) if(.not. myPlasticity(p)) cycle i = i + 1 - NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs + Nconstituents = count(material_phaseAt==p) * discretization_nIPs l = 0 do t = 1,4 do s = 1,param(i)%sum_N_sl @@ -1601,13 +1601,13 @@ end subroutine plastic_nonlocal_results !-------------------------------------------------------------------------------------------------- !> @brief populates the initial dislocation density !-------------------------------------------------------------------------------------------------- -subroutine stateInit(ini,phase,NipcMyPhase,instance) +subroutine stateInit(ini,phase,Nconstituents,instance) type(tInitialParameters) :: & ini integer,intent(in) :: & phase, & - NipcMyPhase, & + Nconstituents, & instance integer :: & e, & @@ -1625,7 +1625,7 @@ subroutine stateInit(ini,phase,NipcMyPhase,instance) totalVolume, & densityBinning, & minimumIpVolume - real(pReal), dimension(NipcMyPhase) :: & + real(pReal), dimension(Nconstituents) :: & volume @@ -1645,13 +1645,13 @@ subroutine stateInit(ini,phase,NipcMyPhase,instance) meanDensity = 0.0_pReal do while(meanDensity < ini%random_rho_u) call random_number(rnd) - phasemember = nint(rnd(1)*real(NipcMyPhase,pReal) + 0.5_pReal) + phasemember = nint(rnd(1)*real(Nconstituents,pReal) + 0.5_pReal) s = nint(rnd(2)*real(sum(ini%N_sl),pReal)*4.0_pReal + 0.5_pReal) meanDensity = meanDensity + densityBinning * volume(phasemember) / totalVolume stt%rhoSglMobile(s,phasemember) = densityBinning enddo else ! homogeneous distribution with noise - do e = 1, NipcMyPhase + do e = 1, Nconstituents do f = 1,size(ini%N_sl,1) from = 1 + sum(ini%N_sl(1:f-1)) upto = sum(ini%N_sl(1:f)) diff --git a/src/constitutive_plastic_phenopowerlaw.f90 b/src/constitutive_plastic_phenopowerlaw.f90 index 50c8e835a..72fa0e9e6 100644 --- a/src/constitutive_plastic_phenopowerlaw.f90 +++ b/src/constitutive_plastic_phenopowerlaw.f90 @@ -70,9 +70,9 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) logical, dimension(:), allocatable :: myPlasticity integer :: & - Ninstance, & + Ninstances, & p, i, & - NipcMyPhase, & + Nconstituents, & sizeState, sizeDotState, & startIndex, endIndex integer, dimension(:), allocatable :: & @@ -91,13 +91,13 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) print'(/,a)', ' <<<+- plastic_phenopowerlaw init -+>>>' myPlasticity = plastic_active('phenopowerlaw') - Ninstance = count(myPlasticity) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myPlasticity) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(dotState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(dotState(Ninstances)) phases => config_material%get('phase') i = 0 @@ -224,20 +224,20 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) !-------------------------------------------------------------------------------------------------- ! allocate state arrays - NipcMyPhase = count(material_phaseAt == p) * discretization_nIPs + Nconstituents = count(material_phaseAt == p) * discretization_nIPs sizeDotState = size(['xi_sl ','gamma_sl']) * prm%sum_N_sl & + size(['xi_tw ','gamma_tw']) * prm%sum_N_tw sizeState = sizeDotState - call constitutive_allocateState(plasticState(p),NipcMyPhase,sizeState,sizeDotState,0) + call constitutive_allocateState(plasticState(p),Nconstituents,sizeState,sizeDotState,0) !-------------------------------------------------------------------------------------------------- ! state aliases and initialization startIndex = 1 endIndex = prm%sum_N_sl stt%xi_slip => plasticState(p)%state (startIndex:endIndex,:) - stt%xi_slip = spread(xi_0_sl, 2, NipcMyPhase) + stt%xi_slip = spread(xi_0_sl, 2, Nconstituents) dot%xi_slip => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) if(any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' @@ -245,7 +245,7 @@ module function plastic_phenopowerlaw_init() result(myPlasticity) startIndex = endIndex + 1 endIndex = endIndex + prm%sum_N_tw stt%xi_twin => plasticState(p)%state (startIndex:endIndex,:) - stt%xi_twin = spread(xi_0_tw, 2, NipcMyPhase) + stt%xi_twin = spread(xi_0_tw, 2, Nconstituents) dot%xi_twin => plasticState(p)%dotState(startIndex:endIndex,:) plasticState(p)%atol(startIndex:endIndex) = pl%get_asFloat('atol_xi',defaultVal=1.0_pReal) if(any(plasticState(p)%atol(startIndex:endIndex) < 0.0_pReal)) extmsg = trim(extmsg)//' atol_xi' diff --git a/src/damage_local.f90 b/src/damage_local.f90 index fefdffce2..e63db90b0 100644 --- a/src/damage_local.f90 +++ b/src/damage_local.f90 @@ -42,7 +42,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine damage_local_init - integer :: Ninstance,NofMyHomog,h + integer :: Ninstances,Nmaterialpoints,h class(tNode), pointer :: & num_generic, & material_homogenization, & @@ -57,8 +57,8 @@ subroutine damage_local_init num%residualStiffness = num_generic%get_asFloat('residualStiffness', defaultVal=1.0e-6_pReal) if (num%residualStiffness < 0.0_pReal) call IO_error(301,ext_msg='residualStiffness') - Ninstance = count(damage_type == DAMAGE_local_ID) - allocate(param(Ninstance)) + Ninstances = count(damage_type == DAMAGE_local_ID) + allocate(param(Ninstances)) material_homogenization => config_material%get('homogenization') do h = 1, material_homogenization%length @@ -73,11 +73,11 @@ subroutine damage_local_init prm%output = homogDamage%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) damageState(h)%sizeState = 1 - allocate(damageState(h)%state0 (1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%subState0(1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%state (1,NofMyHomog), source=damage_initialPhi(h)) + allocate(damageState(h)%state0 (1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%subState0(1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%state (1,Nmaterialpoints), source=damage_initialPhi(h)) nullify(damageMapping(h)%p) damageMapping(h)%p => material_homogenizationMemberAt diff --git a/src/damage_none.f90 b/src/damage_none.f90 index 52100707e..2279bc06b 100644 --- a/src/damage_none.f90 +++ b/src/damage_none.f90 @@ -16,18 +16,18 @@ contains !-------------------------------------------------------------------------------------------------- subroutine damage_none_init - integer :: h,NofMyHomog + integer :: h,Nmaterialpoints print'(/,a)', ' <<<+- damage_none init -+>>>'; flush(6) do h = 1, size(material_name_homogenization) if (damage_type(h) /= DAMAGE_NONE_ID) cycle - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) damageState(h)%sizeState = 0 - allocate(damageState(h)%state0 (0,NofMyHomog)) - allocate(damageState(h)%subState0(0,NofMyHomog)) - allocate(damageState(h)%state (0,NofMyHomog)) + allocate(damageState(h)%state0 (0,Nmaterialpoints)) + allocate(damageState(h)%subState0(0,Nmaterialpoints)) + allocate(damageState(h)%state (0,Nmaterialpoints)) deallocate(damage(h)%p) allocate (damage(h)%p(1), source=damage_initialPhi(h)) diff --git a/src/damage_nonlocal.f90 b/src/damage_nonlocal.f90 index 8b494678e..24a51cf54 100644 --- a/src/damage_nonlocal.f90 +++ b/src/damage_nonlocal.f90 @@ -46,7 +46,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine damage_nonlocal_init - integer :: Ninstance,NofMyHomog,h + integer :: Ninstances,Nmaterialpoints,h class(tNode), pointer :: & num_generic, & material_homogenization, & @@ -60,8 +60,8 @@ subroutine damage_nonlocal_init num_generic => config_numerics%get('generic',defaultVal= emptyDict) num%charLength = num_generic%get_asFloat('charLength',defaultVal=1.0_pReal) - Ninstance = count(damage_type == DAMAGE_nonlocal_ID) - allocate(param(Ninstance)) + Ninstances = count(damage_type == DAMAGE_nonlocal_ID) + allocate(param(Ninstances)) material_homogenization => config_material%get('homogenization') do h = 1, material_homogenization%length @@ -76,11 +76,11 @@ subroutine damage_nonlocal_init prm%output = homogDamage%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) damageState(h)%sizeState = 1 - allocate(damageState(h)%state0 (1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%subState0(1,NofMyHomog), source=damage_initialPhi(h)) - allocate(damageState(h)%state (1,NofMyHomog), source=damage_initialPhi(h)) + allocate(damageState(h)%state0 (1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%subState0(1,Nmaterialpoints), source=damage_initialPhi(h)) + allocate(damageState(h)%state (1,Nmaterialpoints), source=damage_initialPhi(h)) nullify(damageMapping(h)%p) damageMapping(h)%p => material_homogenizationMemberAt diff --git a/src/homogenization_mech_RGC.f90 b/src/homogenization_mech_RGC.f90 index 1a0f7a5d3..5fcea1f8d 100644 --- a/src/homogenization_mech_RGC.f90 +++ b/src/homogenization_mech_RGC.f90 @@ -81,9 +81,9 @@ module subroutine mech_RGC_init(num_homogMech) num_homogMech !< pointer to mechanical homogenization numerics data integer :: & - Ninstance, & + Ninstances, & h, & - NofMyHomog, & + Nmaterialpoints, & sizeState, nIntFaceTot class (tNode), pointer :: & @@ -94,8 +94,8 @@ module subroutine mech_RGC_init(num_homogMech) print'(/,a)', ' <<<+- homogenization_mech_rgc init -+>>>' - Ninstance = count(homogenization_type == HOMOGENIZATION_RGC_ID) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) + Ninstances = count(homogenization_type == HOMOGENIZATION_RGC_ID) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) print*, 'Tjahjanto et al., International Journal of Material Forming 2(1):939–942, 2009' print*, 'https://doi.org/10.1007/s12289-009-0619-1'//IO_EOL @@ -105,10 +105,10 @@ module subroutine mech_RGC_init(num_homogMech) - allocate(param(Ninstance)) - allocate(state(Ninstance)) - allocate(state0(Ninstance)) - allocate(dependentState(Ninstance)) + allocate(param(Ninstances)) + allocate(state(Ninstances)) + allocate(state0(Ninstances)) + allocate(dependentState(Ninstances)) num_RGC => num_homogMech%get('RGC',defaultVal=emptyDict) @@ -173,7 +173,7 @@ module subroutine mech_RGC_init(num_homogMech) prm%D_alpha = homogMech%get_asFloats('D_alpha', requiredSize=3) prm%a_g = homogMech%get_asFloats('a_g', requiredSize=3) - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) nIntFaceTot = 3*( (prm%N_constituents(1)-1)*prm%N_constituents(2)*prm%N_constituents(3) & + prm%N_constituents(1)*(prm%N_constituents(2)-1)*prm%N_constituents(3) & + prm%N_constituents(1)*prm%N_constituents(2)*(prm%N_constituents(3)-1)) @@ -181,24 +181,24 @@ module subroutine mech_RGC_init(num_homogMech) + size(['avg constitutive work ','average penalty energy']) homogState(h)%sizeState = sizeState - allocate(homogState(h)%state0 (sizeState,NofMyHomog), source=0.0_pReal) - allocate(homogState(h)%subState0(sizeState,NofMyHomog), source=0.0_pReal) - allocate(homogState(h)%state (sizeState,NofMyHomog), source=0.0_pReal) + allocate(homogState(h)%state0 (sizeState,Nmaterialpoints), source=0.0_pReal) + allocate(homogState(h)%subState0(sizeState,Nmaterialpoints), source=0.0_pReal) + allocate(homogState(h)%state (sizeState,Nmaterialpoints), source=0.0_pReal) stt%relaxationVector => homogState(h)%state(1:nIntFaceTot,:) st0%relaxationVector => homogState(h)%state0(1:nIntFaceTot,:) stt%work => homogState(h)%state(nIntFaceTot+1,:) stt%penaltyEnergy => homogState(h)%state(nIntFaceTot+2,:) - allocate(dst%volumeDiscrepancy( NofMyHomog), source=0.0_pReal) - allocate(dst%relaxationRate_avg( NofMyHomog), source=0.0_pReal) - allocate(dst%relaxationRate_max( NofMyHomog), source=0.0_pReal) - allocate(dst%mismatch( 3,NofMyHomog), source=0.0_pReal) + allocate(dst%volumeDiscrepancy( Nmaterialpoints), source=0.0_pReal) + allocate(dst%relaxationRate_avg( Nmaterialpoints), source=0.0_pReal) + allocate(dst%relaxationRate_max( Nmaterialpoints), source=0.0_pReal) + allocate(dst%mismatch( 3,Nmaterialpoints), source=0.0_pReal) !-------------------------------------------------------------------------------------------------- ! assigning cluster orientations - dependentState(homogenization_typeInstance(h))%orientation = spread(eu2om(prm%a_g*inRad),3,NofMyHomog) - !dst%orientation = spread(eu2om(prm%a_g*inRad),3,NofMyHomog) ifort version 18.0.1 crashes (for whatever reason) + dependentState(homogenization_typeInstance(h))%orientation = spread(eu2om(prm%a_g*inRad),3,Nmaterialpoints) + !dst%orientation = spread(eu2om(prm%a_g*inRad),3,Nmaterialpoints) ifort version 18.0.1 crashes (for whatever reason) end associate diff --git a/src/homogenization_mech_isostrain.f90 b/src/homogenization_mech_isostrain.f90 index 7f3724ae1..994c1b410 100644 --- a/src/homogenization_mech_isostrain.f90 +++ b/src/homogenization_mech_isostrain.f90 @@ -18,7 +18,7 @@ submodule(homogenization) homogenization_mech_isostrain mapping end type - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -29,9 +29,9 @@ contains module subroutine mech_isostrain_init integer :: & - Ninstance, & + Ninstances, & h, & - NofMyHomog + Nmaterialpoints class(tNode), pointer :: & material_homogenization, & homog, & @@ -39,10 +39,10 @@ module subroutine mech_isostrain_init print'(/,a)', ' <<<+- homogenization_mech_isostrain init -+>>>' - Ninstance = count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) + Ninstances = count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) - allocate(param(Ninstance)) ! one container of parameters per instance + allocate(param(Ninstances)) ! one container of parameters per instance material_homogenization => config_material%get('homogenization') do h = 1, size(homogenization_type) @@ -61,11 +61,11 @@ module subroutine mech_isostrain_init call IO_error(211,ext_msg='sum'//' (mech_isostrain)') end select - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) homogState(h)%sizeState = 0 - allocate(homogState(h)%state0 (0,NofMyHomog)) - allocate(homogState(h)%subState0(0,NofMyHomog)) - allocate(homogState(h)%state (0,NofMyHomog)) + allocate(homogState(h)%state0 (0,Nmaterialpoints)) + allocate(homogState(h)%subState0(0,Nmaterialpoints)) + allocate(homogState(h)%state (0,Nmaterialpoints)) end associate diff --git a/src/homogenization_mech_none.f90 b/src/homogenization_mech_none.f90 index 80949010e..5b12247cd 100644 --- a/src/homogenization_mech_none.f90 +++ b/src/homogenization_mech_none.f90 @@ -14,14 +14,14 @@ contains module subroutine mech_none_init integer :: & - Ninstance, & + Ninstances, & h, & - NofMyHomog + Nmaterialpoints print'(/,a)', ' <<<+- homogenization_mech_none init -+>>>' - Ninstance = count(homogenization_type == HOMOGENIZATION_NONE_ID) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) + Ninstances = count(homogenization_type == HOMOGENIZATION_NONE_ID) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) do h = 1, size(homogenization_type) if(homogenization_type(h) /= HOMOGENIZATION_NONE_ID) cycle @@ -29,11 +29,11 @@ module subroutine mech_none_init if(homogenization_Nconstituents(h) /= 1) & call IO_error(211,ext_msg='N_constituents (mech_none)') - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) homogState(h)%sizeState = 0 - allocate(homogState(h)%state0 (0,NofMyHomog)) - allocate(homogState(h)%subState0(0,NofMyHomog)) - allocate(homogState(h)%state (0,NofMyHomog)) + allocate(homogState(h)%state0 (0,Nmaterialpoints)) + allocate(homogState(h)%subState0(0,Nmaterialpoints)) + allocate(homogState(h)%state (0,Nmaterialpoints)) enddo diff --git a/src/kinematics_cleavage_opening.f90 b/src/kinematics_cleavage_opening.f90 index d52fdbc1c..44bbcb4f4 100644 --- a/src/kinematics_cleavage_opening.f90 +++ b/src/kinematics_cleavage_opening.f90 @@ -20,7 +20,7 @@ submodule(constitutive:constitutive_damage) kinematics_cleavage_opening cleavage_systems end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -35,7 +35,7 @@ module function kinematics_cleavage_opening_init(kinematics_length) result(myKin integer, intent(in) :: kinematics_length logical, dimension(:,:), allocatable :: myKinematics - integer :: Ninstance,p,k + integer :: Ninstances,p,k integer, dimension(:), allocatable :: N_cl !< active number of cleavage systems per family character(len=pStringLen) :: extmsg = '' class(tNode), pointer :: & @@ -48,12 +48,12 @@ module function kinematics_cleavage_opening_init(kinematics_length) result(myKin print'(/,a)', ' <<<+- kinematics_cleavage_opening init -+>>>' myKinematics = kinematics_active('cleavage_opening',kinematics_length) - Ninstance = count(myKinematics) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myKinematics) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(kinematics_cleavage_opening_instance(phases%length), source=0) do p = 1, phases%length diff --git a/src/kinematics_slipplane_opening.f90 b/src/kinematics_slipplane_opening.f90 index e0de5e181..ea8f51427 100644 --- a/src/kinematics_slipplane_opening.f90 +++ b/src/kinematics_slipplane_opening.f90 @@ -22,7 +22,7 @@ submodule(constitutive:constitutive_damage) kinematics_slipplane_opening P_n end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -37,7 +37,7 @@ module function kinematics_slipplane_opening_init(kinematics_length) result(myKi integer, intent(in) :: kinematics_length logical, dimension(:,:), allocatable :: myKinematics - integer :: Ninstance,p,i,k + integer :: Ninstances,p,i,k character(len=pStringLen) :: extmsg = '' integer, dimension(:), allocatable :: N_sl real(pReal), dimension(:,:), allocatable :: d,n,t @@ -51,13 +51,13 @@ module function kinematics_slipplane_opening_init(kinematics_length) result(myKi print'(/,a)', ' <<<+- kinematics_slipplane init -+>>>' myKinematics = kinematics_active('slipplane_opening',kinematics_length) - Ninstance = count(myKinematics) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myKinematics) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') allocate(kinematics_slipplane_opening_instance(phases%length), source=0) - allocate(param(Ninstance)) + allocate(param(Ninstances)) do p = 1, phases%length if(any(myKinematics(:,p))) kinematics_slipplane_opening_instance(p) = count(myKinematics(:,1:p)) diff --git a/src/kinematics_thermal_expansion.f90 b/src/kinematics_thermal_expansion.f90 index 772f5abbf..4afde82e3 100644 --- a/src/kinematics_thermal_expansion.f90 +++ b/src/kinematics_thermal_expansion.f90 @@ -29,7 +29,7 @@ module function kinematics_thermal_expansion_init(kinematics_length) result(myKi integer, intent(in) :: kinematics_length logical, dimension(:,:), allocatable :: myKinematics - integer :: Ninstance,p,i,k + integer :: Ninstances,p,i,k real(pReal), dimension(:), allocatable :: temp class(tNode), pointer :: & phases, & @@ -41,12 +41,12 @@ module function kinematics_thermal_expansion_init(kinematics_length) result(myKi print'(/,a)', ' <<<+- kinematics_thermal_expansion init -+>>>' myKinematics = kinematics_active('thermal_expansion',kinematics_length) - Ninstance = count(myKinematics) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(myKinematics) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(kinematics_thermal_expansion_instance(phases%length), source=0) do p = 1, phases%length diff --git a/src/source_damage_anisoBrittle.f90 b/src/source_damage_anisoBrittle.f90 index faeb8bd87..ca8d6ec2b 100644 --- a/src/source_damage_anisoBrittle.f90 +++ b/src/source_damage_anisoBrittle.f90 @@ -25,7 +25,7 @@ submodule (constitutive:constitutive_damage) source_damage_anisoBrittle output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -45,19 +45,19 @@ module function source_damage_anisoBrittle_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p integer, dimension(:), allocatable :: N_cl character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_anisoBrittle init -+>>>' mySources = source_active('damage_anisoBrittle',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_anisoBrittle_offset (phases%length), source=0) allocate(source_damage_anisoBrittle_instance(phases%length), source=0) @@ -100,8 +100,8 @@ module function source_damage_anisoBrittle_init(source_length) result(mySources) if (any(prm%g_crit < 0.0_pReal)) extmsg = trim(extmsg)//' g_crit' if (any(prm%s_crit < 0.0_pReal)) extmsg = trim(extmsg)//' s_crit' - NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('anisobrittle_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_atol' diff --git a/src/source_damage_anisoDuctile.f90 b/src/source_damage_anisoDuctile.f90 index fd3fa38ed..2fdd46bb4 100644 --- a/src/source_damage_anisoDuctile.f90 +++ b/src/source_damage_anisoDuctile.f90 @@ -19,7 +19,7 @@ submodule(constitutive:constitutive_damage) source_damage_anisoDuctile output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -39,19 +39,19 @@ module function source_damage_anisoDuctile_init(source_length) result(mySources) pl, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p integer, dimension(:), allocatable :: N_sl character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_anisoDuctile init -+>>>' mySources = source_active('damage_anisoDuctile',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_anisoDuctile_offset (phases%length), source=0) allocate(source_damage_anisoDuctile_instance(phases%length), source=0) @@ -84,8 +84,8 @@ module function source_damage_anisoDuctile_init(source_length) result(mySources) if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' if (any(prm%gamma_crit < 0.0_pReal)) extmsg = trim(extmsg)//' gamma_crit' - NipcMyPhase=count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents=count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('anisoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_atol' diff --git a/src/source_damage_isoBrittle.f90 b/src/source_damage_isoBrittle.f90 index c56407c3d..7fcf17ee0 100644 --- a/src/source_damage_isoBrittle.f90 +++ b/src/source_damage_isoBrittle.f90 @@ -17,7 +17,7 @@ submodule(constitutive:constitutive_damage) source_damage_isoBrittle output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -36,18 +36,18 @@ module function source_damage_isoBrittle_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_isoBrittle init -+>>>' mySources = source_active('damage_isoBrittle',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_isoBrittle_offset (phases%length), source=0) allocate(source_damage_isoBrittle_instance(phases%length), source=0) @@ -73,8 +73,8 @@ module function source_damage_isoBrittle_init(source_length) result(mySources) ! sanity checks if (prm%W_crit <= 0.0_pReal) extmsg = trim(extmsg)//' W_crit' - NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,1) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,1) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('isoBrittle_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isobrittle_atol' diff --git a/src/source_damage_isoDuctile.f90 b/src/source_damage_isoDuctile.f90 index 0d068b6e4..1bff20570 100644 --- a/src/source_damage_isoDuctile.f90 +++ b/src/source_damage_isoDuctile.f90 @@ -18,7 +18,7 @@ submodule (constitutive:constitutive_damage) source_damage_isoDuctile output end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -38,18 +38,18 @@ module function source_damage_isoDuctile_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p character(len=pStringLen) :: extmsg = '' print'(/,a)', ' <<<+- source_damage_isoDuctile init -+>>>' mySources = source_active('damage_isoDuctile',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_damage_isoDuctile_offset (phases%length), source=0) allocate(source_damage_isoDuctile_instance(phases%length), source=0) @@ -77,8 +77,8 @@ module function source_damage_isoDuctile_init(source_length) result(mySources) if (prm%q <= 0.0_pReal) extmsg = trim(extmsg)//' q' if (prm%gamma_crit <= 0.0_pReal) extmsg = trim(extmsg)//' gamma_crit' - NipcMyPhase=count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents=count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) sourceState(p)%p(sourceOffset)%atol = src%get_asFloat('isoDuctile_atol',defaultVal=1.0e-3_pReal) if(any(sourceState(p)%p(sourceOffset)%atol < 0.0_pReal)) extmsg = trim(extmsg)//' isoductile_atol' diff --git a/src/source_thermal_dissipation.f90 b/src/source_thermal_dissipation.f90 index 734451a72..f28567aa7 100644 --- a/src/source_thermal_dissipation.f90 +++ b/src/source_thermal_dissipation.f90 @@ -15,7 +15,7 @@ submodule(constitutive:constitutive_thermal) source_thermal_dissipation kappa !< TAYLOR-QUINNEY factor end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -35,17 +35,17 @@ module function source_thermal_dissipation_init(source_length) result(mySources) phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p print'(/,a)', ' <<<+- source_thermal_dissipation init -+>>>' mySources = source_active('thermal_dissipation',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_thermal_dissipation_offset (phases%length), source=0) allocate(source_thermal_dissipation_instance(phases%length), source=0) @@ -61,8 +61,8 @@ module function source_thermal_dissipation_init(source_length) result(mySources) src => sources%get(sourceOffset) prm%kappa = src%get_asFloat('kappa') - NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,0,0,0) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,0,0,0) end associate endif @@ -74,7 +74,7 @@ end function source_thermal_dissipation_init !-------------------------------------------------------------------------------------------------- -!> @brief Ninstances dissipation rate +!> @brief Ninstancess dissipation rate !-------------------------------------------------------------------------------------------------- module subroutine source_thermal_dissipation_getRateAndItsTangent(TDot, dTDot_dT, Tstar, Lp, phase) diff --git a/src/source_thermal_externalheat.f90 b/src/source_thermal_externalheat.f90 index cbc1fa69d..9ba4a051b 100644 --- a/src/source_thermal_externalheat.f90 +++ b/src/source_thermal_externalheat.f90 @@ -19,7 +19,7 @@ submodule(constitutive:constitutive_thermal) source_thermal_externalheat nIntervals end type tParameters - type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstance) + type(tParameters), dimension(:), allocatable :: param !< containers of constitutive parameters (len Ninstances) contains @@ -39,17 +39,17 @@ module function source_thermal_externalheat_init(source_length) result(mySources phase, & sources, & src - integer :: Ninstance,sourceOffset,NipcMyPhase,p + integer :: Ninstances,sourceOffset,Nconstituents,p print'(/,a)', ' <<<+- source_thermal_externalHeat init -+>>>' mySources = source_active('thermal_externalheat',source_length) - Ninstance = count(mySources) - print'(a,i2)', ' # instances: ',Ninstance; flush(IO_STDOUT) - if(Ninstance == 0) return + Ninstances = count(mySources) + print'(a,i2)', ' # instances: ',Ninstances; flush(IO_STDOUT) + if(Ninstances == 0) return phases => config_material%get('phase') - allocate(param(Ninstance)) + allocate(param(Ninstances)) allocate(source_thermal_externalheat_offset (phases%length), source=0) allocate(source_thermal_externalheat_instance(phases%length), source=0) @@ -69,8 +69,8 @@ module function source_thermal_externalheat_init(source_length) result(mySources prm%f_T = src%get_asFloats('f_T',requiredSize = size(prm%t_n)) - NipcMyPhase = count(material_phaseAt==p) * discretization_nIPs - call constitutive_allocateState(sourceState(p)%p(sourceOffset),NipcMyPhase,1,1,0) + Nconstituents = count(material_phaseAt==p) * discretization_nIPs + call constitutive_allocateState(sourceState(p)%p(sourceOffset),Nconstituents,1,1,0) end associate endif diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index 189a7131a..aa807924c 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -40,7 +40,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine thermal_adiabatic_init - integer :: maxNinstance,h,NofMyHomog + integer :: maxNinstances,h,Nmaterialpoints class(tNode), pointer :: & material_homogenization, & homog, & @@ -48,10 +48,10 @@ subroutine thermal_adiabatic_init print'(/,a)', ' <<<+- thermal_adiabatic init -+>>>'; flush(6) - maxNinstance = count(thermal_type == THERMAL_adiabatic_ID) - if (maxNinstance == 0) return + maxNinstances = count(thermal_type == THERMAL_adiabatic_ID) + if (maxNinstances == 0) return - allocate(param(maxNinstance)) + allocate(param(maxNinstances)) material_homogenization => config_material%get('homogenization') do h = 1, size(material_name_homogenization) @@ -67,17 +67,17 @@ subroutine thermal_adiabatic_init prm%output = homogThermal%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog=count(material_homogenizationAt==h) + Nmaterialpoints=count(material_homogenizationAt==h) thermalState(h)%sizeState = 1 - allocate(thermalState(h)%state0 (1,NofMyHomog), source=thermal_initialT(h)) - allocate(thermalState(h)%subState0(1,NofMyHomog), source=thermal_initialT(h)) - allocate(thermalState(h)%state (1,NofMyHomog), source=thermal_initialT(h)) + allocate(thermalState(h)%state0 (1,Nmaterialpoints), source=thermal_initialT(h)) + allocate(thermalState(h)%subState0(1,Nmaterialpoints), source=thermal_initialT(h)) + allocate(thermalState(h)%state (1,Nmaterialpoints), source=thermal_initialT(h)) thermalMapping(h)%p => material_homogenizationMemberAt deallocate(temperature(h)%p) temperature(h)%p => thermalState(h)%state(1,:) deallocate(temperatureRate(h)%p) - allocate (temperatureRate(h)%p(NofMyHomog), source=0.0_pReal) + allocate (temperatureRate(h)%p(Nmaterialpoints), source=0.0_pReal) end associate enddo diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index ead06306a..daa7391a9 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -41,7 +41,7 @@ contains !-------------------------------------------------------------------------------------------------- subroutine thermal_conduction_init - integer :: Ninstance,NofMyHomog,h + integer :: Ninstances,Nmaterialpoints,h class(tNode), pointer :: & material_homogenization, & homog, & @@ -49,8 +49,8 @@ subroutine thermal_conduction_init print'(/,a)', ' <<<+- thermal_conduction init -+>>>'; flush(6) - Ninstance = count(thermal_type == THERMAL_conduction_ID) - allocate(param(Ninstance)) + Ninstances = count(thermal_type == THERMAL_conduction_ID) + allocate(param(Ninstances)) material_homogenization => config_material%get('homogenization') do h = 1, size(material_name_homogenization) @@ -65,17 +65,17 @@ subroutine thermal_conduction_init prm%output = homogThermal%get_asStrings('output',defaultVal=emptyStringArray) #endif - NofMyHomog=count(material_homogenizationAt==h) + Nmaterialpoints=count(material_homogenizationAt==h) thermalState(h)%sizeState = 0 - allocate(thermalState(h)%state0 (0,NofMyHomog)) - allocate(thermalState(h)%subState0(0,NofMyHomog)) - allocate(thermalState(h)%state (0,NofMyHomog)) + allocate(thermalState(h)%state0 (0,Nmaterialpoints)) + allocate(thermalState(h)%subState0(0,Nmaterialpoints)) + allocate(thermalState(h)%state (0,Nmaterialpoints)) thermalMapping(h)%p => material_homogenizationMemberAt deallocate(temperature (h)%p) - allocate (temperature (h)%p(NofMyHomog), source=thermal_initialT(h)) + allocate (temperature (h)%p(Nmaterialpoints), source=thermal_initialT(h)) deallocate(temperatureRate(h)%p) - allocate (temperatureRate(h)%p(NofMyHomog), source=0.0_pReal) + allocate (temperatureRate(h)%p(Nmaterialpoints), source=0.0_pReal) end associate enddo diff --git a/src/thermal_isothermal.f90 b/src/thermal_isothermal.f90 index 703a9aaac..39c8efe91 100644 --- a/src/thermal_isothermal.f90 +++ b/src/thermal_isothermal.f90 @@ -16,18 +16,18 @@ contains !-------------------------------------------------------------------------------------------------- subroutine thermal_isothermal_init - integer :: h,NofMyHomog + integer :: h,Nmaterialpoints print'(/,a)', ' <<<+- thermal_isothermal init -+>>>'; flush(6) do h = 1, size(material_name_homogenization) if (thermal_type(h) /= THERMAL_isothermal_ID) cycle - NofMyHomog = count(material_homogenizationAt == h) + Nmaterialpoints = count(material_homogenizationAt == h) thermalState(h)%sizeState = 0 - allocate(thermalState(h)%state0 (0,NofMyHomog)) - allocate(thermalState(h)%subState0(0,NofMyHomog)) - allocate(thermalState(h)%state (0,NofMyHomog)) + allocate(thermalState(h)%state0 (0,Nmaterialpoints)) + allocate(thermalState(h)%subState0(0,Nmaterialpoints)) + allocate(thermalState(h)%state (0,Nmaterialpoints)) deallocate(temperature (h)%p) allocate (temperature (h)%p(1), source=thermal_initialT(h)) From c3a413e023827f83d969813b79be2e7db8910179 Mon Sep 17 00:00:00 2001 From: Test User Date: Wed, 28 Oct 2020 02:01:59 +0100 Subject: [PATCH 20/36] [skip ci] updated version information after successful test of v3.0.0-alpha-562-gea4904499 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 2e3222238..7cd736d87 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha-556-gf379aecd8 +v3.0.0-alpha-562-gea4904499 From 5f03e8cf8f9a2d22a1157eec311b71286fab4e1a Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 09:31:55 +0100 Subject: [PATCH 21/36] 0-based material indices --- python/damask/_geom.py | 6 +++--- python/tests/test_Geom.py | 1 + 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index f7ac6c437..c9122ca09 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -417,7 +417,7 @@ class Geom: Number of periods per unit cell. Defaults to 1. materials : (int, int), optional Material IDs. Defaults to (1,2). - + Notes ----- The following triply-periodic minimal surfaces are implemented: @@ -687,10 +687,10 @@ class Geom: def renumber(self): - """Renumber sorted material indices to 1,...,N.""" + """Renumber sorted material indices to 0,...,N-1.""" renumbered = np.empty(self.grid,dtype=self.material.dtype) for i, oldID in enumerate(np.unique(self.material)): - renumbered = np.where(self.material == oldID, i+1, renumbered) + renumbered = np.where(self.material == oldID, i, renumbered) return Geom(material = renumbered, size = self.size, diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 368d8ec8f..aa6cf4baf 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -176,6 +176,7 @@ class TestGeom: material = default.material.copy() for m in np.unique(material): material[material==m] = material.max() + np.random.randint(1,30) + default.material -= 1 modified = Geom(material, default.size, default.origin) From b38a498197ad8f20cb099ff1d39bd7538f7ca2ea Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 09:43:20 +0100 Subject: [PATCH 22/36] fast --- python/damask/_geom.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index c9122ca09..c084b8041 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -688,11 +688,9 @@ class Geom: def renumber(self): """Renumber sorted material indices to 0,...,N-1.""" - renumbered = np.empty(self.grid,dtype=self.material.dtype) - for i, oldID in enumerate(np.unique(self.material)): - renumbered = np.where(self.material == oldID, i, renumbered) + _,renumbered = np.unique(self.material,return_inverse=True) - return Geom(material = renumbered, + return Geom(material = renumbered.reshape(self.grid), size = self.size, origin = self.origin, comments = self.comments+[util.execution_stamp('Geom','renumber')], From d72343c865ff7d161c73c3f1f6ad3b31b844c8a3 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 11:41:34 +0100 Subject: [PATCH 23/36] keep order of unique values found in table --- PRIVATE | 2 +- python/damask/_geom.py | 21 +++++++++++++-------- 2 files changed, 14 insertions(+), 9 deletions(-) diff --git a/PRIVATE b/PRIVATE index 1e8c66897..3112a4dbf 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 1e8c66897820468ab46958d995005e2b69204d0e +Subproject commit 3112a4dbfa1e926c07b7f9443161239b8a7e85ca diff --git a/python/damask/_geom.py b/python/damask/_geom.py index c084b8041..d5abe29a7 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -4,6 +4,7 @@ from functools import partial from os import path import numpy as np +import pandas as pd import h5py from scipy import ndimage,spatial @@ -260,15 +261,19 @@ class Geom: Each unique combintation of values results in a material. """ - t = table.sort_by([f'{i}_{coordinates}' for i in range(3,0,-1)]) - - grid,size,origin = grid_filters.cell_coord0_gridSizeOrigin(t.get(coordinates)) + coords = table.sort_by([f'{i}_{coordinates}' for i in range(3,0,-1)]).get(coordinates) + grid,size,origin = grid_filters.cell_coord0_gridSizeOrigin(coords) labels_ = [labels] if isinstance(labels,str) else labels - _,unique_inverse = np.unique(np.hstack([t.get(l) for l in labels_]),return_inverse=True,axis=0) - ma = unique_inverse.reshape(grid,order='F') + 1 + unique,unique_inverse = np.unique(np.hstack([table.get(l) for l in labels_]),return_inverse=True,axis=0) + if len(unique) == grid.prod(): + ma = np.arange(grid.prod()) + else: + ma = np.empty(grid.prod(),'i') + for to_ma,from_ma in enumerate(pd.unique(unique_inverse)): + ma[unique_inverse==from_ma] = to_ma - return Geom(ma,size,origin,util.execution_stamp('Geom','from_table')) + return Geom(ma.reshape(grid,order='F'),size,origin,util.execution_stamp('Geom','from_table')) @staticmethod @@ -782,8 +787,8 @@ class Geom: """ substituted = self.material.copy() - for from_ms,to_ms in zip(from_material,to_material): - substituted[self.material==from_ms] = to_ms + for from_ma,to_ma in zip(from_material,to_material): + substituted[self.material==from_ma] = to_ma return Geom(material = substituted, size = self.size, From 98e0ef3881ca2c13794c741618500ad9f0fd64cd Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 13:38:20 +0100 Subject: [PATCH 24/36] no loops taken from https://stackoverflow.com/questions/3403973 --- python/damask/_geom.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index d5abe29a7..e0e6fcdc6 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -269,9 +269,10 @@ class Geom: if len(unique) == grid.prod(): ma = np.arange(grid.prod()) else: - ma = np.empty(grid.prod(),'i') - for to_ma,from_ma in enumerate(pd.unique(unique_inverse)): - ma[unique_inverse==from_ma] = to_ma + from_ma = pd.unique(unique_inverse) + sort_idx = np.argsort(from_ma) + idx = np.searchsorted(from_ma,unique_inverse,sorter = sort_idx) + ma = np.arange(from_ma.size)[sort_idx][idx] return Geom(ma.reshape(grid,order='F'),size,origin,util.execution_stamp('Geom','from_table')) From 855bf124d397450e988c5e6a1704a5eec5290727 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 14:39:41 +0100 Subject: [PATCH 25/36] faster https://stackoverflow.com/questions/16992713 --- python/damask/_geom.py | 11 +++++++---- python/tests/test_Geom.py | 7 +++++++ 2 files changed, 14 insertions(+), 4 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index e0e6fcdc6..7dfca4cd1 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -787,11 +787,14 @@ class Geom: New material indices. """ - substituted = self.material.copy() - for from_ma,to_ma in zip(from_material,to_material): - substituted[self.material==from_ma] = to_ma + mapper = dict(zip(from_material,to_material)) + + def mp(entry): + return mapper[entry] if entry in mapper else entry + + mp = np.vectorize(mp) - return Geom(material = substituted, + return Geom(material = mp(self.material).reshape(self.grid), size = self.size, origin = self.origin, comments = self.comments+[util.execution_stamp('Geom','substitute')], diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index aa6cf4baf..73d157cfd 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -195,6 +195,13 @@ class TestGeom: modified.substitute(np.arange(default.material.max())+1+offset, np.arange(default.material.max())+1)) + def test_substitute_invariant(self,default): + f = np.unique(default.material.flatten())[:np.random.randint(1,default.material.max())] + t = np.random.permutation(f) + modified = default.substitute(f,t) + assert np.array_equiv(t,f) or (not geom_equal(modified,default)) + assert geom_equal(default, modified.substitute(t,f)) + @pytest.mark.parametrize('axis_angle',[np.array([1,0,0,86.7]), np.array([0,1,0,90.4]), np.array([0,0,1,90]), np.array([1,0,0,175]),np.array([0,-1,0,178]),np.array([0,0,1,180])]) From 76a57af10183861eea045a97cb56783ca6ab941e Mon Sep 17 00:00:00 2001 From: Test User Date: Wed, 28 Oct 2020 19:21:02 +0100 Subject: [PATCH 26/36] [skip ci] updated version information after successful test of v3.0.0-alpha-567-g21f095c9d --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index 7cd736d87..b2c5b72c6 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha-562-gea4904499 +v3.0.0-alpha-567-g21f095c9d From 4ff99a7af4c75cd6c4e7adafb7951ef8134e8a84 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 21:51:20 +0100 Subject: [PATCH 27/36] more logical layout --- python/damask/_geom.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 7dfca4cd1..c76cc8c5a 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -787,14 +787,13 @@ class Geom: New material indices. """ - mapper = dict(zip(from_material,to_material)) - - def mp(entry): + def mp(entry,mapper): return mapper[entry] if entry in mapper else entry mp = np.vectorize(mp) + mapper = dict(zip(from_material,to_material)) - return Geom(material = mp(self.material).reshape(self.grid), + return Geom(material = mp(self.material,mapper).reshape(self.grid), size = self.size, origin = self.origin, comments = self.comments+[util.execution_stamp('Geom','substitute')], From 1b2cd6caf6277ba0b169e799bffadce863b7e6d4 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 28 Oct 2020 21:52:51 +0100 Subject: [PATCH 28/36] documentation was misleading --- python/damask/_configmaterial.py | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index a7ecec108..47b5fa37c 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -229,11 +229,9 @@ class ConfigMaterial(Config): Parameters ---------- constituents : dict - Entries for 'constituents'. The key is the name and the value specifies - the label of the data column in the table + Entries for 'constituents' as key-value pair. **kwargs - Keyword arguments where the key is the name and the value specifies - the label of the data column in the table + Key-value pairs. Examples -------- @@ -264,7 +262,7 @@ class ConfigMaterial(Config): for k,v in kwargs.items(): if hasattr(v,'__len__') and not isinstance(v,str): if len(v) != len(c): - raise ValueError('len mismatch 1') + raise ValueError('Cannot add entries of different length') for i,vv in enumerate(v): c[i][k] = [w.item() for w in vv] if isinstance(vv,np.ndarray) else vv.item() else: @@ -288,7 +286,7 @@ class ConfigMaterial(Config): for k,v in kwargs.items(): if hasattr(v,'__len__') and not isinstance(v,str): if len(v) != N_material: - raise ValueError('len mismatch 2') + raise ValueError('Cannot add entries of different length') for i,vv in enumerate(np.array(v)): m[i][0][k] = [w.item() for w in vv] if isinstance(vv,np.ndarray) else vv.item() else: From 3be0c462a801a628dc556303bc04371e496b4258 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 29 Oct 2020 07:17:41 +0100 Subject: [PATCH 29/36] explicit is better then implicit: user should sort Table according to coordinates to create geometry. This allows to have consistent behavior for from_table in Geom and ConfigMaterial. We always ensure to keep the order --- python/damask/_configmaterial.py | 19 +++++-------------- python/damask/_geom.py | 4 ++-- 2 files changed, 7 insertions(+), 16 deletions(-) diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index 47b5fa37c..c625e6251 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -26,7 +26,7 @@ class ConfigMaterial(Config): @staticmethod - def from_table(table,coordinates=None,constituents={},**kwargs): + def from_table(table,constituents={},**kwargs): """ Load from an ASCII table. @@ -34,10 +34,6 @@ class ConfigMaterial(Config): ---------- table : damask.Table Table that contains material information. - coordinates : str, optional - Label of spatial coordiates. Used for sorting and performing a - sanity check. Default to None, in which case no sorting or checking is - peformed. constituents : dict, optional Entries for 'constituents'. The key is the name and the value specifies the label of the data column in the table @@ -54,7 +50,7 @@ class ConfigMaterial(Config): pos pos pos qu qu qu qu phase homog 0 0 0 0 0.19 0.8 0.24 -0.51 Aluminum SX 1 1 0 0 0.8 0.19 0.24 -0.51 Steel SX - >>> cm.from_table(t,'pos',{'O':'qu','phase':'phase'},homogenization='homog') + >>> cm.from_table(t,{'O':'qu','phase':'phase'},homogenization='homog') material: - constituents: - O: [0.19, 0.8, 0.24, -0.51] @@ -68,17 +64,12 @@ class ConfigMaterial(Config): homogenization: SX """ - if coordinates is not None: - t = table.sort_by([f'{i}_{coordinates}' for i in range(3,0,-1)]) - grid_filters.coord0_check(t.get(coordinates)) - else: - t = table - - constituents_ = {k:t.get(v) for k,v in constituents.items()} - kwargs_ = {k:t.get(v) for k,v in kwargs.items()} + constituents_ = {k:table.get(v) for k,v in constituents.items()} + kwargs_ = {k:table.get(v) for k,v in kwargs.items()} _,idx = np.unique(np.hstack(list({**constituents_,**kwargs_}.values())),return_index=True,axis=0) + idx = np.sort(idx) constituents_ = {k:v[idx].squeeze() for k,v in constituents_.items()} kwargs_ = {k:v[idx].squeeze() for k,v in kwargs_.items()} diff --git a/python/damask/_geom.py b/python/damask/_geom.py index c76cc8c5a..6c5e3907a 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -255,13 +255,13 @@ class Geom: table : damask.Table Table that contains material information. coordinates : str - Label of the column containing the spatial coordinates. + Label of the column containing the vector of spatial coordinates. + Need to be ordered (1./x fast, 3./z slow). labels : str or list of str Label(s) of the columns containing the material definition. Each unique combintation of values results in a material. """ - coords = table.sort_by([f'{i}_{coordinates}' for i in range(3,0,-1)]).get(coordinates) grid,size,origin = grid_filters.cell_coord0_gridSizeOrigin(coords) labels_ = [labels] if isinstance(labels,str) else labels From 999cf53c079c846be8eba6c4d3cbb160db6f2290 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 29 Oct 2020 07:42:41 +0100 Subject: [PATCH 30/36] tests+fixes --- python/damask/_geom.py | 2 +- python/tests/test_ConfigMaterial.py | 13 +++++++++++++ python/tests/test_Geom.py | 13 +++++++++++++ 3 files changed, 27 insertions(+), 1 deletion(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 6c5e3907a..48f022805 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -262,7 +262,7 @@ class Geom: Each unique combintation of values results in a material. """ - grid,size,origin = grid_filters.cell_coord0_gridSizeOrigin(coords) + grid,size,origin = grid_filters.cell_coord0_gridSizeOrigin(table.get(coordinates)) labels_ = [labels] if isinstance(labels,str) else labels unique,unique_inverse = np.unique(np.hstack([table.get(l) for l in labels_]),return_inverse=True,axis=0) diff --git a/python/tests/test_ConfigMaterial.py b/python/tests/test_ConfigMaterial.py index 4863a8ac4..d0bcd2a56 100644 --- a/python/tests/test_ConfigMaterial.py +++ b/python/tests/test_ConfigMaterial.py @@ -1,8 +1,10 @@ import os import pytest +import numpy as np from damask import ConfigMaterial +from damask import Table @pytest.fixture def reference_dir(reference_dir_base): @@ -74,3 +76,14 @@ class TestConfigMaterial: material_config = ConfigMaterial.load(reference_dir/'material.yaml') new = material_config.material_rename_homogenization({'Taylor':'isostrain'}) assert not new.is_complete + + def test_from_table(self): + N = np.random.randint(3,10) + a = np.vstack((np.hstack((np.arange(N),np.arange(N)[::-1])),np.ones(N*2),np.zeros(N*2),np.ones(N*2))).T + t = Table(a,{'varying':2,'constant':2}) + c = ConfigMaterial.from_table(t,constituents={'a':'varying','b':'1_constant'},c='2_constant') + assert len(c['material']) == N + for i,m in enumerate(c['material']): + c = m['constituents'][0] + assert m['c'] == 1 and c['b'] == 0 and c['a'] == [i,1] + diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 73d157cfd..dd9c28a3e 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -3,8 +3,10 @@ import numpy as np from damask import VTK from damask import Geom +from damask import Table from damask import Rotation from damask import util +from damask import grid_filters def geom_equal(a,b): @@ -371,3 +373,14 @@ class TestGeom: grid = np.ones(3,dtype=int)*64 geom = Geom.from_minimal_surface(grid,np.ones(3),surface,threshold) assert np.isclose(np.count_nonzero(geom.material==1)/np.prod(geom.grid),.5,rtol=1e-3) + + + def test_from_table(self): + grid = np.random.randint(60,100,3) + size = np.ones(3)+np.random.rand(3) + coords = grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') + x = np.ones(grid.prod()).reshape(-1,1) + y = np.hstack([np.arange(grid[1])]*grid[0]*grid[2]).reshape(-1,1) + t = Table(np.hstack((coords,x,y)),{'coords':3,'x':1,'y':1}) + g = Geom.from_table(t,'coords',['x','y']) + assert g.N_materials == g.grid[2] From 2736bbaba7a6931a01d9981a6b206977e5ec51e4 Mon Sep 17 00:00:00 2001 From: Test User Date: Thu, 29 Oct 2020 15:49:28 +0100 Subject: [PATCH 31/36] [skip ci] updated version information after successful test of v3.0.0-alpha-579-g4b25097a9 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index b2c5b72c6..eaf9a5b74 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha-567-g21f095c9d +v3.0.0-alpha-579-g4b25097a9 From 4b700d367e1ebb6e07453e167348028fcb840d29 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 29 Oct 2020 17:52:16 +0100 Subject: [PATCH 32/36] improved test --- python/damask/_configmaterial.py | 1 - python/tests/test_Geom.py | 11 ++++++----- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index c625e6251..8b4f1c5ea 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -2,7 +2,6 @@ import copy import numpy as np -from . import grid_filters from . import Config from . import Lattice from . import Rotation diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index dd9c28a3e..03bd92e13 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -379,8 +379,9 @@ class TestGeom: grid = np.random.randint(60,100,3) size = np.ones(3)+np.random.rand(3) coords = grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') - x = np.ones(grid.prod()).reshape(-1,1) - y = np.hstack([np.arange(grid[1])]*grid[0]*grid[2]).reshape(-1,1) - t = Table(np.hstack((coords,x,y)),{'coords':3,'x':1,'y':1}) - g = Geom.from_table(t,'coords',['x','y']) - assert g.N_materials == g.grid[2] + z=np.ones(grid.prod()) + z[int(z.size/2):]=0 + t = Table(np.column_stack((coords,z)),{'coords':3,'z':1}) + g = Geom.from_table(t,'coords',['1_coords','z']) + assert g.N_materials == g.grid[0]*2 and \ + (g.material[:,:,-1]-g.material[:,:,0]-grid[:1].prod() == 0).all() From 0e499eedf367d64e2fba6c4dc13532c31ce579b8 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 29 Oct 2020 18:35:22 +0100 Subject: [PATCH 33/36] correct rounding/clearer logic --- python/tests/test_Geom.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 03bd92e13..0d393a4c4 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -380,8 +380,7 @@ class TestGeom: size = np.ones(3)+np.random.rand(3) coords = grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') z=np.ones(grid.prod()) - z[int(z.size/2):]=0 + z[grid[:2].prod()*int(grid[2]/2):]=0 t = Table(np.column_stack((coords,z)),{'coords':3,'z':1}) g = Geom.from_table(t,'coords',['1_coords','z']) - assert g.N_materials == g.grid[0]*2 and \ - (g.material[:,:,-1]-g.material[:,:,0]-grid[:1].prod() == 0).all() + assert g.N_materials == g.grid[0]*2 and (g.material[:,:,-1]-g.material[:,:,0] == grid[0]).all() From b2289e80b22c1fa9eceee5362752a885e0259e74 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 29 Oct 2020 20:09:13 +0100 Subject: [PATCH 34/36] simplified --- python/damask/_configmaterial.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index 8b4f1c5ea..0f8d4efa1 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -226,10 +226,11 @@ class ConfigMaterial(Config): Examples -------- >>> import damask - >>> m = damask.ConfigMaterial() >>> O = damask.Rotation.from_random(3).as_quaternion() >>> phase = ['Aluminum','Steel','Aluminum'] - >>> m.material_add(constituents={'phase':phase,'O':O},homogenization='SX') + >>> m = damask.ConfigMaterial().material_add(constituents={'phase':phase,'O':O}, + ... homogenization='SX') + >>> m material: - constituents: - O: [0.577764, -0.146299, -0.617669, 0.513010] @@ -259,8 +260,7 @@ class ConfigMaterial(Config): for i in range(len(c)): c[i][k] = v dup = copy.deepcopy(self) - if 'material' not in dup: dup['material'] = [] - dup['material'] +=c + dup['material'] = dup['material'] + c if 'material' in dup else c return dup From 10ba53721dec264e9d94d215df6a4b13294e124d Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 29 Oct 2020 20:47:15 +0100 Subject: [PATCH 35/36] not needed anymore --- examples/.gitignore | 3 --- 1 file changed, 3 deletions(-) diff --git a/examples/.gitignore b/examples/.gitignore index 3cb6b9820..93d78295b 100644 --- a/examples/.gitignore +++ b/examples/.gitignore @@ -3,6 +3,3 @@ *.xdmf *.sta *.vt* -*.geom -*.seeds -postProc From 0af81e1af7e808d46695f4431ead4d0692390359 Mon Sep 17 00:00:00 2001 From: Test User Date: Fri, 30 Oct 2020 05:57:50 +0100 Subject: [PATCH 36/36] [skip ci] updated version information after successful test of v3.0.0-alpha-594-g46e5023f8 --- VERSION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/VERSION b/VERSION index eaf9a5b74..002deaea7 100644 --- a/VERSION +++ b/VERSION @@ -1 +1 @@ -v3.0.0-alpha-579-g4b25097a9 +v3.0.0-alpha-594-g46e5023f8