From 27d861decc0fc6d5d7d96a427029833e5b05eef1 Mon Sep 17 00:00:00 2001
From: Martin	Diehl <m.diehl@mpie.de>
Date: Thu, 3 Jul 2014 13:17:29 +0000
Subject: [PATCH] fixed calculation of post results size introduced with new
 state

---
 code/Makefile                 |  4 +-
 code/constitutive.f90         | 97 ++++++++++++-----------------------
 code/constitutive_damage.f90  | 56 +++++++++-----------
 code/constitutive_thermal.f90 | 56 +++++++++-----------
 code/crystallite.f90          |  1 -
 5 files changed, 82 insertions(+), 132 deletions(-)

diff --git a/code/Makefile b/code/Makefile
index 903502940..59b04f278 100644
--- a/code/Makefile
+++ b/code/Makefile
@@ -29,7 +29,7 @@ endif
 COMPILERNAME ?= $(F90)
 
 INCLUDE_DIRS :=-I../lib
-LIBRARIES :=-lfftw3
+LIBRARIES :=-llapack -lfftw3
 LIB_DIRS  :=-L$(FFTW_ROOT)/lib
 RUN_PATH  :=-Wl,-rpath,$(FFTW_ROOT)/lib
 
@@ -372,7 +372,7 @@ endif
 DAMASK_spectral.exe: DAMASK_spectral_driver.o 
 	$(PREFIX) $(COMPILERNAME) $(OPENMP_FLAG_$(F90)) $(LINK_OPTIONS_$(F90)) $(STANDARD_CHECK_$(F90))  $(OPTIMIZATION_$(MAXOPTI)_$(F90)) \
   -o DAMASK_spectral.exe DAMASK_spectral_driver.o \
-  $(COMPILED_FILES) $(LIBRARIES) $(LIB_DIRS) $(RUN_PATH) $(SUFFIX)
+  $(COMPILED_FILES) $(LIB_DIRS) $(LIBRARIES) $(LIB_DIRS) $(RUN_PATH) $(SUFFIX)
 
 DAMASK_spectral_driver.o: DAMASK_spectral_driver.f90 DAMASK_spectral_solverBasic.o $(PETSC_FILES)
 	$(PREFIX) $(COMPILERNAME) $(COMPILE_MAXOPTI) -c DAMASK_spectral_driver.f90 $(SUFFIX)
diff --git a/code/constitutive.f90 b/code/constitutive.f90
index f668a680d..31c19bd6b 100644
--- a/code/constitutive.f90
+++ b/code/constitutive.f90
@@ -11,8 +11,6 @@ module constitutive
  
  implicit none
  private
- integer(pInt), public, dimension(:,:,:), allocatable :: &
-   constitutive_sizePostResults                                                                      !< size of postResults array per grain
  integer(pInt), public, protected :: &
    constitutive_maxSizePostResults, &
    constitutive_maxSizeDotState
@@ -104,17 +102,9 @@ subroutine constitutive_init
  implicit none
  integer(pInt), parameter :: FILEUNIT = 200_pInt
  integer(pInt) :: &
-  g, &                                                                                              !< grain number
-  i, &                                                                                              !< integration point number
-  e, &                                                                                              !< element number
-  cMax, &                                                                                           !< maximum number of grains
-  iMax, &                                                                                           !< maximum number of integration points
-  eMax, &                                                                                           !< maximum number of elements
+  e, &                                                                                        !< maximum number of elements
   phase, &
-  s, &
-  p, &
-  instance,&
-  myNgrains
+  instance
 
  integer(pInt), dimension(:,:), pointer :: thisSize
  character(len=64), dimension(:,:), pointer :: thisOutput
@@ -185,54 +175,41 @@ subroutine constitutive_init
  enddo
  close(FILEUNIT)
  
-!--------------------------------------------------------------------------------------------------
-! allocation of states
- cMax = homogenization_maxNgrains
- iMax = mesh_maxNips
- eMax = mesh_NcpElems
- allocate(constitutive_sizePostResults(cMax,iMax,eMax), source=0_pInt) 
- ElemLoop:do e = 1_pInt,mesh_NcpElems                                                               ! loop over elements
-   myNgrains = homogenization_Ngrains(mesh_element(3,e)) 
-   IPloop:do i = 1_pInt,FE_Nips(FE_geomtype(mesh_element(2,e)))                                     ! loop over IPs
-     GrainLoop:do g = 1_pInt,myNgrains                                                              ! loop over grains
-       select case(phase_elasticity(material_phase(g,i,e)))                                            
-         case default                                                                               ! so far no output for elasticity
-       end select
-       phase = material_phase(g,i,e)
-       instance = phase_plasticityInstance(phase)
-       select case(phase_plasticity(material_phase(g,i,e)))
-         case (PLASTICITY_NONE_ID)
-           constitutive_sizePostResults(g,i,e) =    0_pInt
-         case (PLASTICITY_J2_ID) 
-           constitutive_sizePostResults(g,i,e) =    constitutive_j2_sizePostResults(instance)
-         case (PLASTICITY_PHENOPOWERLAW_ID)
-           constitutive_sizePostResults(g,i,e) =    constitutive_phenopowerlaw_sizePostResults(instance)
-         case (PLASTICITY_DISLOTWIN_ID)
-           constitutive_sizePostResults(g,i,e) =    constitutive_dislotwin_sizePostResults(instance)
-         case (PLASTICITY_TITANMOD_ID)
-           constitutive_sizePostResults(g,i,e) =   constitutive_titanmod_sizePostResults(instance)
-         case (PLASTICITY_NONLOCAL_ID)
-           nonlocalConstitutionPresent = .true.
-           plasticState(mappingConstitutive(2,g,i,e))%nonlocal = .true.
-           if(myNgrains/=1_pInt) call IO_error(252_pInt, e,i,g)
-           constitutive_sizePostResults(g,i,e) =    constitutive_nonlocal_sizePostResults(instance)
-       end select
-     enddo GrainLoop
-   enddo IPloop
- enddo ElemLoop
+
+ PhaseLoop:do phase = 1_pInt,material_Nphase                                                              ! loop over phases
+   instance = phase_plasticityInstance(phase)
+   select case(phase_plasticity(phase))
+     case (PLASTICITY_NONE_ID)
+       plasticState(phase)%sizePostResults = constitutive_none_sizePostResults(instance)
+     case (PLASTICITY_J2_ID)
+       plasticState(phase)%sizePostResults = constitutive_j2_sizePostResults(instance)
+     case (PLASTICITY_PHENOPOWERLAW_ID)
+       plasticState(phase)%sizePostResults = constitutive_none_sizePostResults(instance)
+     case (PLASTICITY_DISLOTWIN_ID)
+       plasticState(phase)%sizePostResults = constitutive_dislotwin_sizePostResults(instance)
+     case (PLASTICITY_TITANMOD_ID)
+       plasticState(phase)%sizePostResults = constitutive_titanmod_sizePostResults(instance)
+     case (PLASTICITY_NONLOCAL_ID)
+       nonlocalConstitutionPresent = .true.
+       plasticState(phase)%nonlocal = .true.
+       plasticState(phase)%sizePostResults = constitutive_nonlocal_sizePostResults(instance)
+   end select   
+
+ enddo PhaseLoop
 
  if (nonlocalConstitutionPresent) &
    call constitutive_nonlocal_stateInit()
 
- do e = 1_pInt,mesh_NcpElems                                                                        ! loop over elements
-   myNgrains = homogenization_Ngrains(mesh_element(3,e)) 
-   forall(i = 1_pInt:FE_Nips(FE_geomtype(mesh_element(2,e))), g = 1_pInt:myNgrains)
-     plasticState(mappingConstitutive(2,g,i,e))%partionedState0(:,mappingConstitutive(1,g,i,e)) = &
-      plasticState(mappingConstitutive(2,g,i,e))%State0(:,mappingConstitutive(1,g,i,e))    ! need to be defined for first call of constitutive_microstructure in crystallite_init
-     plasticState(mappingConstitutive(2,g,i,e))%State(:,mappingConstitutive(1,g,i,e)) =          &
-      plasticState(mappingConstitutive(2,g,i,e))%State0(:,mappingConstitutive(1,g,i,e))    ! need to be defined for first call of constitutive_microstructure in crystallite_init
-   endforall
- enddo
+
+ constitutive_maxSizeDotState = 0_pInt
+ constitutive_maxSizePostResults = 0_pInt
+
+ PhaseLoop2:do phase = 1_pInt,material_Nphase
+   plasticState(phase)%partionedState0 = plasticState(phase)%State0
+   plasticState(phase)%State = plasticState(phase)%State0
+   constitutive_maxSizeDotState = max(constitutive_maxSizeDotState, plasticState(phase)%sizeDotState)
+   constitutive_maxSizePostResults = max(constitutive_maxSizePostResults, plasticState(phase)%sizePostResults)
+ enddo PhaseLoop2
 
 #ifdef HDF
  call  HDF5_mappingConstitutive(mappingConstitutive)
@@ -276,12 +253,6 @@ subroutine constitutive_init
  flush(6)
 #endif
 
- constitutive_maxSizePostResults = 0_pInt
- constitutive_maxSizeDotState = 0_pInt
- do p = 1, size(plasticState)
-  constitutive_maxSizeDotState = max(constitutive_maxSizeDotState, plasticState(p)%sizeDotState)
-  constitutive_maxSizePostResults = max(constitutive_maxSizePostResults, plasticState(p)%sizePostResults)
- enddo
 
 end subroutine constitutive_init
 
@@ -733,7 +704,7 @@ function constitutive_postResults(Tstar_v, FeArray, temperature, ipc, ip, el)
    ipc, &                                                                                           !< grain number
    ip, &                                                                                            !< integration point number
    el                                                                                               !< element number
- real(pReal), dimension(constitutive_sizePostResults(ipc,ip,el)) :: &
+ real(pReal), dimension(plasticState(mappingConstitutive(1,ipc,ip,el))%sizePostResults) :: &
    constitutive_postResults
  real(pReal),  intent(in) :: &
    temperature
diff --git a/code/constitutive_damage.f90 b/code/constitutive_damage.f90
index 837c62e1c..1395f726f 100644
--- a/code/constitutive_damage.f90
+++ b/code/constitutive_damage.f90
@@ -63,17 +63,9 @@ use damage_gradient
  implicit none
  integer(pInt), parameter :: FILEUNIT = 200_pInt
  integer(pInt) :: &
-  g, &                                                                                              !< grain number
-  i, &                                                                                              !< integration point number
-  e, &                                                                                              !< element number
-  cMax, &                                                                                           !< maximum number of grains
-  iMax, &                                                                                           !< maximum number of integration points
-  eMax, &                                                                                           !< maximum number of elements
-  phase, &
-  s, &
-  p, &
-  instance,&
-  myNgrains
+  e, &                                                                                              !< grain number
+  ph, &
+  instance
 
  integer(pInt), dimension(:,:), pointer :: thisSize
  logical :: knownDamage
@@ -96,10 +88,10 @@ use damage_gradient
 !--------------------------------------------------------------------------------------------------
 ! write description file for constitutive phase output
  call IO_write_jobFile(FILEUNIT,'outputDamage') 
- do phase = 1_pInt,material_Nphase
-   instance = phase_damageInstance(phase)                                                           ! which instance of a plasticity is present phase
+ do ph = 1_pInt,material_Nphase
+   instance = phase_damageInstance(ph)                                                           ! which instance of a plasticity is present phase
    knownDamage = .true.
-   select case(phase_damage(phase))                                                                 ! split per constititution
+   select case(phase_damage(ph))                                                                 ! split per constititution
      case (DAMAGE_none_ID)
        outputName = DAMAGE_NONE_label
        thisOutput => null()
@@ -111,11 +103,11 @@ use damage_gradient
      case default
        knownDamage = .false.
    end select   
-   write(FILEUNIT,'(/,a,/)') '['//trim(phase_name(phase))//']'
+   write(FILEUNIT,'(/,a,/)') '['//trim(phase_name(ph))//']'
    if (knownDamage) then
      write(FILEUNIT,'(a)') '(damage)'//char(9)//trim(outputName)
-     if (phase_damage(phase) /= DAMAGE_none_ID) then
-       do e = 1_pInt,phase_Noutput(phase)
+     if (phase_damage(ph) /= DAMAGE_none_ID) then
+       do e = 1_pInt,phase_Noutput(ph)
          write(FILEUNIT,'(a,i4)') trim(thisOutput(e,instance))//char(9),thisSize(e,instance)
        enddo
      endif
@@ -125,24 +117,22 @@ use damage_gradient
  
 !--------------------------------------------------------------------------------------------------
 ! allocation of states
- PhaseLoop:do phase = 1_pInt,material_Nphase                                                              ! loop over phases
-   instance = phase_damageInstance(phase)
-   select case(phase_damage(phase))
-     case (DAMAGE_none_ID) 
-       damageState(phase)%sizePostResults = damage_none_sizePostResults(instance)
-
-     case (DAMAGE_gradient_ID) 
-       damageState(phase)%sizePostResults = damage_gradient_sizePostResults(instance)
-       
-   end select
- enddo PhaseLoop
- 
  constitutive_damage_maxSizePostResults = 0_pInt
  constitutive_damage_maxSizeDotState = 0_pInt
- do p = 1, size(damageState)
-  constitutive_damage_maxSizeDotState = max(constitutive_damage_maxSizeDotState, damageState(p)%sizeDotState)
-  constitutive_damage_maxSizePostResults = max(constitutive_damage_maxSizePostResults, damageState(p)%sizePostResults)
- enddo
+ PhaseLoop:do ph = 1_pInt,material_Nphase                                                              ! loop over phases
+   instance = phase_damageInstance(ph)
+   select case(phase_damage(ph))
+     case (DAMAGE_none_ID) 
+       damageState(ph)%sizePostResults = damage_none_sizePostResults(instance)
+
+     case (DAMAGE_gradient_ID) 
+       damageState(ph)%sizePostResults = damage_gradient_sizePostResults(instance)
+       
+   end select
+  constitutive_damage_maxSizeDotState = max(constitutive_damage_maxSizeDotState, damageState(ph)%sizeDotState)
+  constitutive_damage_maxSizePostResults = max(constitutive_damage_maxSizePostResults, damageState(ph)%sizePostResults)
+ enddo PhaseLoop
+
 end subroutine constitutive_damage_init
 
 
diff --git a/code/constitutive_thermal.f90 b/code/constitutive_thermal.f90
index 3b80e7b76..37f8f3ab0 100644
--- a/code/constitutive_thermal.f90
+++ b/code/constitutive_thermal.f90
@@ -63,17 +63,9 @@ subroutine constitutive_thermal_init
  implicit none
  integer(pInt), parameter :: FILEUNIT = 200_pInt
  integer(pInt) :: &
-  g, &                                                                                              !< grain number
-  i, &                                                                                              !< integration point number
-  e, &                                                                                              !< element number
-  cMax, &                                                                                           !< maximum number of grains
-  iMax, &                                                                                           !< maximum number of integration points
-  eMax, &                                                                                           !< maximum number of elements
-  phase, &
-  s, &
-  p, &
-  instance,&
-  myNgrains
+  e, &                                                                                              !< grain number
+  ph, &                                                                                                !< phase
+  instance
 
  integer(pInt), dimension(:,:), pointer :: thisSize
  logical :: knownThermal
@@ -96,10 +88,10 @@ subroutine constitutive_thermal_init
 !--------------------------------------------------------------------------------------------------
 ! write description file for constitutive phase output
  call IO_write_jobFile(FILEUNIT,'outputThermal') 
- do phase = 1_pInt,material_Nphase
-   instance = phase_thermalInstance(phase)                                                           ! which instance is present phase
+ do ph = 1_pInt,material_Nphase
+   instance = phase_thermalInstance(ph)                                                              ! which instance is present phase
    knownThermal = .true.
-   select case(phase_thermal(phase))                                                                 ! split per constititution
+   select case(phase_thermal(ph))                                                                 ! split per constititution
      case (THERMAL_none_ID)
        outputName = THERMAL_NONE_label
        thisOutput => null()
@@ -111,11 +103,11 @@ subroutine constitutive_thermal_init
      case default
        knownThermal = .false.
    end select   
-   write(FILEUNIT,'(/,a,/)') '['//trim(phase_name(phase))//']'
+   write(FILEUNIT,'(/,a,/)') '['//trim(phase_name(ph))//']'
    if (knownThermal) then
      write(FILEUNIT,'(a)') '(thermal)'//char(9)//trim(outputName)
-     if (phase_thermal(phase) /= THERMAL_none_ID) then
-       do e = 1_pInt,phase_Noutput(phase)
+     if (phase_thermal(ph) /= THERMAL_none_ID) then
+       do e = 1_pInt,phase_Noutput(ph)
          write(FILEUNIT,'(a,i4)') trim(thisOutput(e,instance))//char(9),thisSize(e,instance)
        enddo
      endif
@@ -125,24 +117,22 @@ subroutine constitutive_thermal_init
  
 !--------------------------------------------------------------------------------------------------
 ! allocation of states
- PhaseLoop:do phase = 1_pInt,material_Nphase                                                              ! loop over phases
-   instance = phase_thermalInstance(phase)
-   select case(phase_thermal(phase))
-     case (THERMAL_none_ID) 
-       thermalState(phase)%sizePostResults = thermal_none_sizePostResults(instance)
-
-     case (THERMAL_conduction_ID) 
-       thermalState(phase)%sizePostResults = thermal_conduction_sizePostResults(instance)
-       
-   end select
- enddo PhaseLoop
- 
  constitutive_thermal_maxSizePostResults = 0_pInt
  constitutive_thermal_maxSizeDotState = 0_pInt
- do p = 1, size(thermalState)
-  constitutive_thermal_maxSizeDotState = max(constitutive_thermal_maxSizeDotState, thermalState(p)%sizeDotState)
-  constitutive_thermal_maxSizePostResults = max(constitutive_thermal_maxSizePostResults, thermalState(p)%sizePostResults)
- enddo
+ PhaseLoop:do ph = 1_pInt,material_Nphase                                                              ! loop over phases
+   instance = phase_thermalInstance(ph)
+   select case(phase_thermal(ph))
+     case (THERMAL_none_ID) 
+       thermalState(ph)%sizePostResults = thermal_none_sizePostResults(instance)
+
+     case (THERMAL_conduction_ID) 
+       thermalState(ph)%sizePostResults = thermal_conduction_sizePostResults(instance)
+       
+   end select
+  constitutive_thermal_maxSizeDotState = max(constitutive_thermal_maxSizeDotState, thermalState(ph)%sizeDotState)
+  constitutive_thermal_maxSizePostResults = max(constitutive_thermal_maxSizePostResults, thermalState(ph)%sizePostResults)
+ enddo PhaseLoop
+
 end subroutine constitutive_thermal_init
 
 
diff --git a/code/crystallite.f90 b/code/crystallite.f90
index 070bf5ed6..7f8558b54 100644
--- a/code/crystallite.f90
+++ b/code/crystallite.f90
@@ -3948,7 +3948,6 @@ function crystallite_postResults(ipc, ip, el)
    material_texture, &
    homogenization_Ngrains
  use constitutive, only: &
-   constitutive_sizePostResults, &
    constitutive_postResults, &
    constitutive_homogenizedC
  use constitutive_damage, only: &