reworked distribution of microstructure constituents. now each material point represents "as good as possible" the fractional content of constituents.

removed error for volume fractions not equalling 1.

implemented capability to rotate the texture given in material.config.
This commit is contained in:
Philip Eisenlohr 2013-05-29 17:23:49 +00:00
parent 32a16f9745
commit 84c43741a6
1 changed files with 189 additions and 109 deletions

View File

@ -116,7 +116,7 @@ module material
material_volume, & !< volume of each grain,IP,element
texture_Gauss, & !< data of each Gauss component
texture_Fiber, & !< data of each Fiber component
texture_Rotation !< rotation of each texture
texture_rotation !< rotation of each texture
logical, dimension(:), allocatable, private :: &
homogenization_active
@ -196,22 +196,22 @@ subroutine material_init
if (minval(microstructure_texture(1:microstructure_Nconstituents(m),m)) < 1_pInt .or. &
maxval(microstructure_texture(1:microstructure_Nconstituents(m),m)) > material_Ntexture) &
call IO_error(152_pInt,m)
if (abs(sum(microstructure_fraction(:,m)) - 1.0_pReal) >= 1.0e-10_pReal) then
if (iand(myDebug,debug_levelExtensive) /= 0_pInt) then
write(6,'(a)') ' sum of microstructure fraction = ',sum(microstructure_fraction(:,m))
endif
call IO_error(153_pInt,m)
endif
! if (abs(sum(microstructure_fraction(:,m)) - 1.0_pReal) >= 1.0e-6_pReal) then ! have ppm precision in fractions
! if (iand(myDebug,debug_levelExtensive) /= 0_pInt) then
! write(6,'(a,1x,f12.9)') ' sum of microstructure fraction = ',sum(microstructure_fraction(:,m))
! endif
! call IO_error(153_pInt,m)
! endif
enddo
debugOut: if (iand(myDebug,debug_levelExtensive) /= 0_pInt) then
write(6,'(/,a,/)') ' MATERIAL configuration'
write(6,'(a32,1x,a16,1x,a6)') 'homogenization ','type ','grains'
do h = 1_pInt,material_Nhomogenization
write(6,'(1x,a32,1x,a16,1x,i4)') homogenization_name(h),homogenization_type(h),homogenization_Ngrains(h)
write(6,'(1x,a32,1x,a16,1x,i6)') homogenization_name(h),homogenization_type(h),homogenization_Ngrains(h)
enddo
write(6,'(/,a32,19x,a11,1x,a12,1x,a13)') 'microstructure','crystallite','constituents','homogeneous'
write(6,'(/,a14,18x,1x,a11,1x,a12,1x,a13)') 'microstructure','crystallite','constituents','homogeneous'
do m = 1_pInt,material_Nmicrostructure
write(6,'(a32,4x,i4,8x,i4,8x,l1)') microstructure_name(m), &
write(6,'(1x,a32,1x,i11,1x,i12,1x,l13)') microstructure_name(m), &
microstructure_crystallite(m), &
microstructure_Nconstituents(m), &
microstructure_elemhomo(m)
@ -603,9 +603,9 @@ subroutine material_parseTexture(myFile,myPart)
texture_maxNfiber = maxval(texture_Nfiber)
allocate(texture_Gauss (5,texture_maxNgauss,Nsections)); texture_Gauss = 0.0_pReal
allocate(texture_Fiber (6,texture_maxNfiber,Nsections)); texture_Fiber = 0.0_pReal
allocate(texture_Rotation(3,3,Nsections));
do j =1_pInt, Nsections
texture_Rotation(1:3,1:3,j) = math_I3
allocate(texture_rotation(3,3,Nsections));
do j = 1_pInt, Nsections
texture_rotation(1:3,1:3,j) = math_I3
enddo
rewind(myFile)
@ -636,23 +636,23 @@ subroutine material_parseTexture(myFile,myPart)
textureType: select case(tag)
case ('rotation') textureType
do j = 1_pInt, 3_pInt
do j = 1_pInt, 3_pInt ! look for "x", "y", and "z" entries
tag = IO_lc(IO_stringValue(line,positions,j+1_pInt))
select case (tag)
case('x', '+x')
texture_Rotation(1_pInt:3_pInt,j,section) = (/1.0_pReal, 0.0_pReal, 0.0_pReal/) ! original axis is now +x-axis
texture_rotation(j,1:3,section) = (/ 1.0_pReal, 0.0_pReal, 0.0_pReal/) ! original axis is now +x-axis
case('-x')
texture_Rotation(1_pInt:3_pInt,j,section) = (/-1.0_pReal, 0.0_pReal, 0.0_pReal/) ! original axis is now -x-axis
texture_rotation(j,1:3,section) = (/-1.0_pReal, 0.0_pReal, 0.0_pReal/) ! original axis is now -x-axis
case('y', '+y')
texture_Rotation(1_pInt:3_pInt,j,section) = (/0.0_pReal, 1.0_pReal, 0.0_pReal/) ! original axis is now +y-axis
texture_rotation(j,1:3,section) = (/ 0.0_pReal, 1.0_pReal, 0.0_pReal/) ! original axis is now +y-axis
case('-y')
texture_Rotation(1_pInt:3_pInt,j,section) = (/0.0_pReal, -1.0_pReal, 0.0_pReal/) ! original axis is now -y-axis
texture_rotation(j,1:3,section) = (/ 0.0_pReal,-1.0_pReal, 0.0_pReal/) ! original axis is now -y-axis
case('z', '+z')
texture_Rotation(1_pInt:3_pInt,j,section) = (/0.0_pReal, 0.0_pReal, 1.0_pReal/) ! original axis is now +z-axis
texture_rotation(j,1:3,section) = (/ 0.0_pReal, 0.0_pReal, 1.0_pReal/) ! original axis is now +z-axis
case('-z')
texture_Rotation(1_pInt:3_pInt,j,section) = (/0.0_pReal, 0.0_pReal, -1.0_pReal/) ! original axis is now -z-axis
texture_rotation(j,1:3,section) = (/ 0.0_pReal, 0.0_pReal,-1.0_pReal/) ! original axis is now -z-axis
case default
call IO_error(156_pInt)
call IO_error(157_pInt,section)
end select
enddo
@ -735,6 +735,10 @@ subroutine material_parseTexture(myFile,myPart)
!--------------------------------------------------------------------------------------------------
subroutine material_populateGrains
use math, only: &
math_RtoEuler, &
math_EulerToR, &
math_mul33x33, &
math_range, &
math_sampleRandomOri, &
math_sampleGaussOri, &
math_sampleFiberOri, &
@ -758,16 +762,20 @@ subroutine material_populateGrains
implicit none
integer(pInt), dimension (:,:), allocatable :: Ngrains
integer(pInt), dimension (microstructure_maxNconstituents) &
:: NgrainsOfConstituent
integer(pInt), dimension (microstructure_maxNconstituents) :: &
NgrainsOfConstituent, &
currentGrainOfConstituent, &
randomOrder
real(pReal), dimension (microstructure_maxNconstituents) :: &
rndArray
real(pReal), dimension (:), allocatable :: volumeOfGrain
real(pReal), dimension (:,:), allocatable :: orientationOfGrain
real(pReal), dimension (3) :: orientation
real(pReal), dimension (3,3) :: symOrientation
integer(pInt), dimension (:), allocatable :: phaseOfGrain, textureOfGrain
integer(pInt) :: t,e,i,g,j,m,homog,micro,sgn,hme, myDebug
integer(pInt) :: phaseID,textureID,dGrains,myNgrains,myNorientations, &
grain,constituentGrain,symExtension
integer(pInt) :: t,e,i,g,j,m,c,r,homog,micro,sgn,hme, myDebug, &
phaseID,textureID,dGrains,myNgrains,myNorientations,myNconstituents, &
grain,constituentGrain,ipGrain,symExtension, ip
real(pReal) :: extreme,rnd
integer(pInt), dimension (:,:), allocatable :: Nelems ! counts number of elements in homog, micro array
type(p_intvec), dimension (:,:), allocatable :: elemsOfHomogMicro ! lists element number in homog, micro array
@ -803,23 +811,23 @@ subroutine material_populateGrains
! identify maximum grain count per IP (from element) and find grains per homog/micro pair
Nelems = 0_pInt ! reuse as counter
elementLooping: do e = 1_pInt,mesh_NcpElems
t = FE_geomtype(mesh_element(2,e))
t = FE_geomtype(mesh_element(2,e))
homog = mesh_element(3,e)
micro = mesh_element(4,e)
if (homog < 1_pInt .or. homog > material_Nhomogenization) & ! out of bounds
call IO_error(154_pInt,e,0_pInt,0_pInt)
if (micro < 1_pInt .or. micro > material_Nmicrostructure) & ! out of bounds
call IO_error(155_pInt,e,0_pInt,0_pInt)
if (microstructure_elemhomo(micro)) then
dGrains = homogenization_Ngrains(homog)
if (microstructure_elemhomo(micro)) then ! how many grains are needed at this element?
dGrains = homogenization_Ngrains(homog) ! only one set of Ngrains (other IPs are plain copies)
else
dGrains = homogenization_Ngrains(homog) * FE_Nips(t)
dGrains = homogenization_Ngrains(homog) * FE_Nips(t) ! each IP has Ngrains
endif
Ngrains(homog,micro) = Ngrains(homog,micro) + dGrains
Nelems(homog,micro) = Nelems(homog,micro) + 1_pInt
Ngrains(homog,micro) = Ngrains(homog,micro) + dGrains ! total grain count
Nelems(homog,micro) = Nelems(homog,micro) + 1_pInt ! total element count
elemsOfHomogMicro(homog,micro)%p(Nelems(homog,micro)) = e ! remember elements active in this homog/micro pair
enddo elementLooping
allocate(volumeOfGrain(maxval(Ngrains))) ! reserve memory for maximum case
allocate(phaseOfGrain(maxval(Ngrains))) ! reserve memory for maximum case
allocate(textureOfGrain(maxval(Ngrains))) ! reserve memory for maximum case
@ -836,41 +844,58 @@ subroutine material_populateGrains
do micro = 1_pInt,material_Nmicrostructure ! all pairs of homog and micro
if (Ngrains(homog,micro) > 0_pInt) then ! an active pair of homog and micro
myNgrains = Ngrains(homog,micro) ! assign short name for total number of grains to populate
myNconstituents = microstructure_Nconstituents(micro) ! assign short name for number of constituents
if (iand(myDebug,debug_levelBasic) /= 0_pInt) then
!$OMP CRITICAL (write2out)
write(6,'(/,a32,1x,a32,1x,i6)') homogenization_name(homog),microstructure_name(micro),myNgrains
!$OMP END CRITICAL (write2out)
endif
!--------------------------------------------------------------------------------------------------
! calculate volume of each grain
volumeOfGrain = 0.0_pReal
grain = 0_pInt
do hme = 1_pInt, Nelems(homog,micro)
e = elemsOfHomogMicro(homog,micro)%p(hme) ! my combination of homog and micro, only perform calculations for elements with homog, micro combinations which is indexed in cpElemsindex
e = elemsOfHomogMicro(homog,micro)%p(hme) ! my combination of homog and micro, only perform calculations for elements with homog, micro combinations which is indexed in cpElemsindex
t = FE_geomtype(mesh_element(2,e))
if (microstructure_elemhomo(micro)) then ! homogeneous distribution of grains over each element's IPs
volumeOfGrain(grain+1_pInt:grain+dGrains) = sum(mesh_ipVolume(1:FE_Nips(t),e))/&
real(dGrains,pReal)
grain = grain + dGrains ! wind forward by NgrainsPerIP
real(dGrains,pReal) ! each grain combines size of all IPs in that element
grain = grain + dGrains ! wind forward by Ngrains@IP
else
forall (i = 1_pInt:FE_Nips(t)) & ! loop over IPs
volumeOfGrain(grain+(i-1)*dGrains+1_pInt:grain+i*dGrains) = &
mesh_ipVolume(i,e)/dGrains ! assign IPvolume/Ngrains to all grains of IP
grain = grain + FE_Nips(t) * dGrains ! wind forward by Nips*NgrainsPerIP
mesh_ipVolume(i,e)/dGrains ! assign IPvolume/Ngrains@IP to all grains of IP
grain = grain + FE_Nips(t) * dGrains ! wind forward by Nips*Ngrains@IP
endif
enddo
if (grain /= myNgrains) &
call IO_error(0,e = homog,i = micro,ext_msg = 'inconsistent grain count after volume calc')
!--------------------------------------------------------------------------------------------------
! divide myNgrains as best over constituents
NgrainsOfConstituent = 0_pInt
forall (i = 1_pInt:microstructure_Nconstituents(micro)) &
!
! example: three constituents with fractions of 0.25, 0.25, and 0.5 distributed over 20 (microstructure) grains
!
! ***** ***** **********
! NgrainsOfConstituent: 5, 5, 10
! counters:
! |-----> grain (if constituent == 2)
! |--> constituentGrain (of constituent 2)
!
NgrainsOfConstituent = 0_pInt ! reset counter of grains per constituent
forall (i = 1_pInt:myNconstituents) &
NgrainsOfConstituent(i) = nint(microstructure_fraction(i,micro) * myNgrains, pInt) ! do rounding integer conversion
do while (sum(NgrainsOfConstituent) /= myNgrains) ! total grain count over constituents wrong?
sgn = sign(1_pInt, myNgrains - sum(NgrainsOfConstituent)) ! direction of required change
extreme = 0.0_pReal
t = 0_pInt
do i = 1_pInt,microstructure_Nconstituents(micro) ! find largest deviator
do i = 1_pInt,myNconstituents ! find largest deviator
if (real(sgn,pReal)*log(NgrainsOfConstituent(i)/myNgrains/microstructure_fraction(i,micro)) > extreme) then
extreme = real(sgn,pReal)*log(NgrainsOfConstituent(i)/myNgrains/microstructure_fraction(i,micro))
t = i
@ -879,116 +904,171 @@ subroutine material_populateGrains
NgrainsOfConstituent(t) = NgrainsOfConstituent(t) + sgn ! change that by one
enddo
!--------------------------------------------------------------------------------------------------
! assign phase and texture info
phaseOfGrain = 0_pInt
textureOfGrain = 0_pInt
orientationOfGrain = 0.0_pReal
grain = 0_pInt ! reset microstructure grain index
constituents: do i = 1_pInt,microstructure_Nconstituents(micro) ! loop over constituents
texture: do i = 1_pInt,myNconstituents ! loop over constituents
grain = sum(NgrainsOfConstituent(1_pInt:i-1_pInt)) ! set microstructure grain index of current constituent
! "grain" points to start of this constituent's grain population
constituentGrain = 0_pInt ! constituent grain index
phaseID = microstructure_phase(i,micro)
textureID = microstructure_texture(i,micro)
phaseOfGrain(grain+1_pInt:grain+NgrainsOfConstituent(i)) = phaseID ! assign resp. phase
phaseOfGrain (grain+1_pInt:grain+NgrainsOfConstituent(i)) = phaseID ! assign resp. phase
textureOfGrain(grain+1_pInt:grain+NgrainsOfConstituent(i)) = textureID ! assign resp. texture
myNorientations = ceiling(real(NgrainsOfConstituent(i),pReal)/&
real(texture_symmetry(textureID),pReal),pInt) ! max number of unique orientations (excl. symmetry)
constituentGrain = 0_pInt ! constituent grain index
!--------------------------------------------------------------------------------------------------
! dealing with texture components
! ...has texture components
if (texture_ODFfile(textureID) == '') then
do t = 1_pInt,texture_Ngauss(textureID) ! loop over Gauss components
gauss: do t = 1_pInt,texture_Ngauss(textureID) ! loop over Gauss components
do g = 1_pInt,int(myNorientations*texture_Gauss(5,t,textureID),pInt) ! loop over required grain count
orientationOfGrain(:,grain+constituentGrain+g) = &
math_sampleGaussOri(texture_Gauss(1:3,t,textureID),&
texture_Gauss( 4,t,textureID))
enddo
constituentGrain = constituentGrain + int(myNorientations*texture_Gauss(5,t,textureID))
enddo
constituentGrain = &
constituentGrain + int(myNorientations*texture_Gauss(5,t,textureID)) ! advance counter for grains of current constituent
enddo gauss
do t = 1_pInt,texture_Nfiber(textureID) ! loop over fiber components
fiber: do t = 1_pInt,texture_Nfiber(textureID) ! loop over fiber components
do g = 1_pInt,int(myNorientations*texture_Fiber(6,t,textureID),pInt) ! loop over required grain count
orientationOfGrain(:,grain+constituentGrain+g) = &
math_sampleFiberOri(texture_Fiber(1:2,t,textureID),&
texture_Fiber(3:4,t,textureID),&
texture_Fiber( 5,t,textureID))
enddo
constituentGrain = constituentGrain + int(myNorientations*texture_fiber(6,t,textureID),pInt)
enddo
do j = constituentGrain+1_pInt,myNorientations ! fill remainder with random
orientationOfGrain(:,grain+j) = math_sampleRandomOri()
enddo
else
constituentGrain = &
constituentGrain + int(myNorientations*texture_fiber(6,t,textureID),pInt) ! advance counter for grains of current constituent
enddo fiber
random: do constituentGrain = constituentGrain+1_pInt,myNorientations ! fill remainder with random
orientationOfGrain(:,grain+constituentGrain) = math_sampleRandomOri()
enddo random
!--------------------------------------------------------------------------------------------------
! hybrid IA
orientationOfGrain(:,grain+1:grain+myNorientations) = IO_hybridIA(myNorientations,texture_ODFfile(textureID))
if (all(orientationOfGrain(:,grain+1) == -1.0_pReal)) call IO_error(156_pInt)
constituentGrain = constituentGrain + myNorientations
! ...has hybrid IA
else
orientationOfGrain(1:3,grain+1_pInt:grain+myNorientations) = &
IO_hybridIA(myNorientations,texture_ODFfile(textureID))
if (all(orientationOfGrain(1:3,grain+1_pInt) == -1.0_pReal)) call IO_error(156_pInt)
endif
!--------------------------------------------------------------------------------------------------
! ...texture rotation
do j = 1_pInt,myNorientations ! loop over each "real" orientation
orientationOfGrain(1:3,grain+j) = math_RtoEuler( & ! translate back to Euler angles
math_mul33x33( & ! pre-multiply
math_EulertoR(orientationOfGrain(1:3,grain+j)), & ! face-value orientation
texture_rotation(1:3,1:3,textureID) & ! rotation matrix and
) &
)
enddo
!--------------------------------------------------------------------------------------------------
! ...sample symmetry
symExtension = texture_symmetry(textureID) - 1_pInt
if (symExtension > 0_pInt) then ! sample symmetry
constituentGrain = NgrainsOfConstituent(i)-myNorientations ! calc remainder of array
if (symExtension > 0_pInt) then ! sample symmetry (number of additional equivalent orientations)
constituentGrain = myNorientations ! start right after "real" orientations
do j = 1_pInt,myNorientations ! loop over each "real" orientation
symOrientation = math_symmetricEulers(texture_symmetry(textureID),orientationOfGrain(:,j)) ! get symmetric equivalents
e = min(symExtension,constituentGrain) ! are we at end of constituent grain array?
symOrientation = math_symmetricEulers(texture_symmetry(textureID), &
orientationOfGrain(1:3,grain+j)) ! get symmetric equivalents
e = min(symExtension,NgrainsOfConstituent(i)-constituentGrain) ! do not overshoot end of constituent grain array
if (e > 0_pInt) then
orientationOfGrain(:,grain+myNorientations+1+(j-1_pInt)*symExtension:&
grain+myNorientations+e+(j-1_pInt)*symExtension) = &
symOrientation(:,1:e)
constituentGrain = constituentGrain - e ! remainder shrinks by e
orientationOfGrain(1:3,grain+constituentGrain+1: &
grain+constituentGrain+e) = &
symOrientation(1:3,1:e)
constituentGrain = constituentGrain + e ! remainder shrinks by e
endif
enddo
endif
grain = grain + NgrainsOfConstituent(i) ! advance microstructure grain index
enddo constituents
!--------------------------------------------------------------------------------------------------
! unless element homogeneous, reshuffle grains
if (.not. microstructure_elemhomo(micro)) then
do i=1_pInt,myNgrains-1_pInt ! walk thru grains
! shuffle grains within current constituent
do j = 1_pInt,NgrainsOfConstituent(i)-1_pInt ! walk thru grains of current constituent
call random_number(rnd)
t = nint(rnd*(myNgrains-i)+i+0.5_pReal,pInt) ! select a grain in remaining list
m = phaseOfGrain(t) ! exchange current with random
phaseOfGrain(t) = phaseOfGrain(i)
phaseOfGrain(i) = m
m = textureOfGrain(t) ! exchange current with random
textureOfGrain(t) = textureOfGrain(i)
textureOfGrain(i) = m
orientation = orientationOfGrain(:,t)
orientationOfGrain(:,t) = orientationOfGrain(:,i)
orientationOfGrain(:,i) = orientation
t = nint(rnd*(NgrainsOfConstituent(i)-j)+j+0.5_pReal,pInt) ! select a grain in remaining list
m = phaseOfGrain(grain+t) ! exchange current with random
phaseOfGrain(grain+t) = phaseOfGrain(grain+j)
phaseOfGrain(grain+j) = m
m = textureOfGrain(grain+t) ! exchange current with random
textureOfGrain(grain+t) = textureOfGrain(grain+j)
textureOfGrain(grain+j) = m
orientation = orientationOfGrain(1:3,grain+t) ! exchange current with random
orientationOfGrain(1:3,grain+t) = orientationOfGrain(1:3,grain+j)
orientationOfGrain(1:3,grain+j) = orientation
enddo
endif
enddo texture
!--------------------------------------------------------------------------------------------------
! calc fraction after weighing with volumePerGrain, exchange in MC steps to improve result
grain = 0_pInt
!>@ToDo calc fraction after weighing with volumePerGrain, exchange in MC steps to improve result (humbug at the moment)
!--------------------------------------------------------------------------------------------------
! distribute grains of all constituents as accurately as possible to given constituent fractions
ip = 0_pInt
currentGrainOfConstituent = 0_pInt
do hme = 1_pInt, Nelems(homog,micro)
e = elemsOfHomogMicro(homog,micro)%p(hme) ! only perform calculations for elements with homog, micro combinations which is indexed in cpElemsindex
t = FE_geomtype(mesh_element(2,e))
if (microstructure_elemhomo(micro)) then ! homogeneous distribution of grains over each element's IPs
forall (i = 1_pInt:FE_Nips(t), g = 1_pInt:dGrains) ! loop over IPs and grains
material_volume(g,i,e) = volumeOfGrain(grain+g)
material_phase(g,i,e) = phaseOfGrain(grain+g)
material_texture(g,i,e) = textureOfGrain(grain+g)
material_EulerAngles(:,g,i,e) = orientationOfGrain(:,grain+g)
end forall
FEsolving_execIP(2,e) = 1_pInt ! restrict calculation to first IP only, since all other results are to be copied from this
grain = grain + dGrains ! wind forward by NgrainsPerIP
m = 1_pInt ! process only first IP
else
forall (i = 1_pInt:FE_Nips(t), g = 1_pInt:dGrains) ! loop over IPs and grains
material_volume(g,i,e) = volumeOfGrain(grain+(i-1_pInt)*dGrains+g)
material_phase(g,i,e) = phaseOfGrain(grain+(i-1_pInt)*dGrains+g)
material_texture(g,i,e) = textureOfGrain(grain+(i-1_pInt)*dGrains+g)
material_EulerAngles(:,g,i,e) = orientationOfGrain(:,grain+(i-1_pInt)*dGrains+g)
end forall
grain = grain + FE_Nips(t) * dGrains ! wind forward by Nips*NgrainsPerIP
m = FE_Nips(t) ! process all IPs
endif
do i = 1_pInt, m ! loop over necessary IPs
ip = ip + 1_pInt ! keep track of total ip count
ipGrain = 0_pInt ! count number of grains assigned at this IP
randomOrder = math_range(microstructure_maxNconstituents) ! start out with ordered sequence of constituents
call random_number(rndArray) ! as many rnd numbers as (max) constituents
do j = 1_pInt, myNconstituents - 1_pInt ! loop over constituents ...
r = nint(rndArray(j)*(myNconstituents-j)+j+0.5_pReal,pInt) ! ... select one in remaining list
c = randomOrder(r) ! ... call it "c"
randomOrder(r) = randomOrder(j) ! ... and exchange with present position in constituent list
grain = sum(NgrainsOfConstituent(1:c-1_pInt)) ! figure out actual starting index in overall/consecutive grain population
do g = 1_pInt, min(dGrains-ipGrain, & ! leftover number of grains at this IP
max(0_pInt, & ! no negative values
nint(real(ip * dGrains * NgrainsOfConstituent(c)) / & ! fraction of grains scaled to this constituent...
real(myNgrains),pInt) - & ! ...minus those already distributed
currentGrainOfConstituent(c)))
ipGrain = ipGrain + 1_pInt ! advance IP grain counter
currentGrainOfConstituent(c) = currentGrainOfConstituent(c) + 1_pInt ! advance index of grain population for constituent c
material_volume(ipGrain,i,e) = volumeOfGrain(grain+currentGrainOfConstituent(c)) ! assign properties
material_phase(ipGrain,i,e) = phaseOfGrain(grain+currentGrainOfConstituent(c))
material_texture(ipGrain,i,e) = textureOfGrain(grain+currentGrainOfConstituent(c))
material_EulerAngles(1:3,ipGrain,i,e) = orientationOfGrain(1:3,grain+currentGrainOfConstituent(c))
enddo; enddo
c = randomOrder(microstructure_Nconstituents(micro)) ! look up constituent remaining after random shuffling
grain = sum(NgrainsOfConstituent(1:c-1_pInt)) ! figure out actual starting index in overall/consecutive grain population
do ipGrain = ipGrain + 1_pInt, dGrains ! ensure last constituent fills up to dGrains
currentGrainOfConstituent(c) = currentGrainOfConstituent(c) + 1_pInt
material_volume(ipGrain,i,e) = volumeOfGrain(grain+currentGrainOfConstituent(c))
material_phase(ipGrain,i,e) = phaseOfGrain(grain+currentGrainOfConstituent(c))
material_texture(ipGrain,i,e) = textureOfGrain(grain+currentGrainOfConstituent(c))
material_EulerAngles(1:3,ipGrain,i,e) = orientationOfGrain(1:3,grain+currentGrainOfConstituent(c))
enddo
enddo
do i = i, FE_Nips(t) ! loop over IPs to (possibly) distribute copies from first IP
material_volume (1_pInt:dGrains,i,e) = material_volume (1_pInt:dGrains,1,e)
material_phase (1_pInt:dGrains,i,e) = material_phase (1_pInt:dGrains,1,e)
material_texture(1_pInt:dGrains,i,e) = material_texture(1_pInt:dGrains,1,e)
material_EulerAngles(1:3,1_pInt:dGrains,i,e) = material_EulerAngles(1:3,1_pInt:dGrains,1,e)
enddo
enddo
endif ! active homog,micro pair
enddo
@ -999,7 +1079,7 @@ subroutine material_populateGrains
deallocate(textureOfGrain)
deallocate(orientationOfGrain)
deallocate(Nelems)
!> ToDo - causing segmentation fault: needs looking into
!>@ToDo - causing segmentation fault: needs looking into
!do homog = 1,material_Nhomogenization
! do micro = 1,material_Nmicrostructure
! if (Nelems(homog,micro) > 0_pInt) deallocate(elemsOfHomogMicro(homog,micro)%p)