diff --git a/code/IO.f90 b/code/IO.f90 index c03789705..50978f120 100644 --- a/code/IO.f90 +++ b/code/IO.f90 @@ -1056,6 +1056,10 @@ endfunction msg = 'Convergence not reached' case (610) msg = 'Stress loop not converged' + + case (666) + msg = 'Memory leak detected' + case (700) msg = 'Singular matrix in stress iteration' case default diff --git a/code/constitutive.f90 b/code/constitutive.f90 index 6f866ed1e..83df3b9b8 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -55,7 +55,7 @@ subroutine constitutive_init() integer(pInt), parameter :: fileunit = 200 integer(pInt) e,i,g,p,myInstance,myNgrains integer(pInt), dimension(:,:), pointer :: thisSize - character(64), dimension(:,:), pointer :: thisOutput + character(len=64), dimension(:,:), pointer :: thisOutput logical knownConstitution if(.not. IO_open_file(fileunit,material_configFile)) call IO_error (100) ! corrupt config file @@ -65,6 +65,7 @@ subroutine constitutive_init() call constitutive_dislotwin_init(fileunit) call constitutive_nonlocal_init(fileunit) + close(fileunit) ! write description file for constitutive phase output diff --git a/code/constitutive_j2.f90 b/code/constitutive_j2.f90 index cf820c250..27044440f 100644 --- a/code/constitutive_j2.f90 +++ b/code/constitutive_j2.f90 @@ -83,6 +83,9 @@ subroutine constitutive_j2_init(file) maxNinstance = count(phase_constitution == constitutive_j2_label) if (maxNinstance == 0) return + write(6,'(a16,x,i5)') '# instances:',maxNinstance + write(6,*) + allocate(constitutive_j2_sizeDotState(maxNinstance)) ; constitutive_j2_sizeDotState = 0_pInt allocate(constitutive_j2_sizeState(maxNinstance)) ; constitutive_j2_sizeState = 0_pInt allocate(constitutive_j2_sizePostResults(maxNinstance)); constitutive_j2_sizePostResults = 0_pInt diff --git a/code/constitutive_phenopowerlaw.f90 b/code/constitutive_phenopowerlaw.f90 index f4cf30f14..3d07d3da3 100644 --- a/code/constitutive_phenopowerlaw.f90 +++ b/code/constitutive_phenopowerlaw.f90 @@ -145,11 +145,13 @@ subroutine constitutive_phenopowerlaw_init(file) write(6,'(a20,a20,a12)') '<<<+- constitutive_',constitutive_phenopowerlaw_label,' init -+>>>' write(6,*) '$Id$' write(6,*) - maxNinstance = count(phase_constitution == constitutive_phenopowerlaw_label) if (maxNinstance == 0) return - + + write(6,'(a16,x,i5)') '# instances:',maxNinstance + write(6,*) + allocate(constitutive_phenopowerlaw_sizeDotState(maxNinstance)) ; constitutive_phenopowerlaw_sizeDotState = 0_pInt allocate(constitutive_phenopowerlaw_sizeState(maxNinstance)) ; constitutive_phenopowerlaw_sizeState = 0_pInt allocate(constitutive_phenopowerlaw_sizePostResults(maxNinstance)); constitutive_phenopowerlaw_sizePostResults = 0_pInt @@ -217,6 +219,7 @@ subroutine constitutive_phenopowerlaw_init(file) rewind(file) line = '' section = 0 + do while (IO_lc(IO_getTag(line,'<','>')) /= 'phase') ! wind forward to read(file,'(a1024)',END=100) line enddo @@ -308,6 +311,7 @@ subroutine constitutive_phenopowerlaw_init(file) enddo 100 do i = 1,maxNinstance + constitutive_phenopowerlaw_structure(i) = lattice_initializeStructure(constitutive_phenopowerlaw_structureName(i), & ! get structure constitutive_phenopowerlaw_CoverA(i)) constitutive_phenopowerlaw_Nslip(:,i) = min(lattice_NslipSystem(:,constitutive_phenopowerlaw_structure(i)),& ! limit active slip systems per family to max @@ -352,33 +356,34 @@ subroutine constitutive_phenopowerlaw_init(file) constitutive_phenopowerlaw_hardeningMatrix_sliptwin = 0.0_pReal constitutive_phenopowerlaw_hardeningMatrix_twinslip = 0.0_pReal constitutive_phenopowerlaw_hardeningMatrix_twintwin = 0.0_pReal + do i = 1,maxNinstance do j = 1,maxval(phase_Noutput) - select case(constitutive_phenopowerlaw_output(j,i)) - case('resistance_slip', & - 'shearrate_slip', & - 'resolvedstress_slip' & - ) - mySize = constitutive_phenopowerlaw_totalNslip(i) - case('resistance_twin', & - 'shearrate_twin', & - 'resolvedstress_twin' & - ) - mySize = constitutive_phenopowerlaw_totalNtwin(i) - case('totalshear', & - 'totalvolfrac' & - ) - mySize = 1_pInt - case default - mySize = 0_pInt - end select + select case(constitutive_phenopowerlaw_output(j,i)) + case('resistance_slip', & + 'shearrate_slip', & + 'resolvedstress_slip' & + ) + mySize = constitutive_phenopowerlaw_totalNslip(i) + case('resistance_twin', & + 'shearrate_twin', & + 'resolvedstress_twin' & + ) + mySize = constitutive_phenopowerlaw_totalNtwin(i) + case('totalshear', & + 'totalvolfrac' & + ) + mySize = 1_pInt + case default + mySize = 0_pInt + end select - if (mySize > 0_pInt) then ! any meaningful output found - constitutive_phenopowerlaw_sizePostResult(j,i) = mySize - constitutive_phenopowerlaw_sizePostResults(i) = & - constitutive_phenopowerlaw_sizePostResults(i) + mySize - endif + if (mySize > 0_pInt) then ! any meaningful output found + constitutive_phenopowerlaw_sizePostResult(j,i) = mySize + constitutive_phenopowerlaw_sizePostResults(i) = & + constitutive_phenopowerlaw_sizePostResults(i) + mySize + endif enddo constitutive_phenopowerlaw_sizeDotState(i) = constitutive_phenopowerlaw_totalNslip(i)+ & @@ -387,27 +392,27 @@ subroutine constitutive_phenopowerlaw_init(file) constitutive_phenopowerlaw_totalNtwin(i)+ 2 ! s_slip, s_twin, sum(gamma), sum(f) select case (constitutive_phenopowerlaw_structure(i)) ! assign elasticity tensor - case(1:2) ! cubic(s) - forall(k=1:3) - forall(j=1:3) & - constitutive_phenopowerlaw_Cslip_66(k,j,i) = constitutive_phenopowerlaw_C12(i) - constitutive_phenopowerlaw_Cslip_66(k,k,i) = constitutive_phenopowerlaw_C11(i) - constitutive_phenopowerlaw_Cslip_66(k+3,k+3,i) = constitutive_phenopowerlaw_C44(i) - end forall - case(3) ! hex - constitutive_phenopowerlaw_Cslip_66(1,1,i) = constitutive_phenopowerlaw_C11(i) - constitutive_phenopowerlaw_Cslip_66(2,2,i) = constitutive_phenopowerlaw_C11(i) - constitutive_phenopowerlaw_Cslip_66(3,3,i) = constitutive_phenopowerlaw_C33(i) - constitutive_phenopowerlaw_Cslip_66(1,2,i) = constitutive_phenopowerlaw_C12(i) - constitutive_phenopowerlaw_Cslip_66(2,1,i) = constitutive_phenopowerlaw_C12(i) - constitutive_phenopowerlaw_Cslip_66(1,3,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(3,1,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(2,3,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(3,2,i) = constitutive_phenopowerlaw_C13(i) - constitutive_phenopowerlaw_Cslip_66(4,4,i) = constitutive_phenopowerlaw_C44(i) - constitutive_phenopowerlaw_Cslip_66(5,5,i) = constitutive_phenopowerlaw_C44(i) - constitutive_phenopowerlaw_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_phenopowerlaw_C11(i)- & - constitutive_phenopowerlaw_C12(i)) + case(1:2) ! cubic(s) + forall(k=1:3) + forall(j=1:3) & + constitutive_phenopowerlaw_Cslip_66(k,j,i) = constitutive_phenopowerlaw_C12(i) + constitutive_phenopowerlaw_Cslip_66(k,k,i) = constitutive_phenopowerlaw_C11(i) + constitutive_phenopowerlaw_Cslip_66(k+3,k+3,i) = constitutive_phenopowerlaw_C44(i) + end forall + case(3) ! hex + constitutive_phenopowerlaw_Cslip_66(1,1,i) = constitutive_phenopowerlaw_C11(i) + constitutive_phenopowerlaw_Cslip_66(2,2,i) = constitutive_phenopowerlaw_C11(i) + constitutive_phenopowerlaw_Cslip_66(3,3,i) = constitutive_phenopowerlaw_C33(i) + constitutive_phenopowerlaw_Cslip_66(1,2,i) = constitutive_phenopowerlaw_C12(i) + constitutive_phenopowerlaw_Cslip_66(2,1,i) = constitutive_phenopowerlaw_C12(i) + constitutive_phenopowerlaw_Cslip_66(1,3,i) = constitutive_phenopowerlaw_C13(i) + constitutive_phenopowerlaw_Cslip_66(3,1,i) = constitutive_phenopowerlaw_C13(i) + constitutive_phenopowerlaw_Cslip_66(2,3,i) = constitutive_phenopowerlaw_C13(i) + constitutive_phenopowerlaw_Cslip_66(3,2,i) = constitutive_phenopowerlaw_C13(i) + constitutive_phenopowerlaw_Cslip_66(4,4,i) = constitutive_phenopowerlaw_C44(i) + constitutive_phenopowerlaw_Cslip_66(5,5,i) = constitutive_phenopowerlaw_C44(i) + constitutive_phenopowerlaw_Cslip_66(6,6,i) = 0.5_pReal*(constitutive_phenopowerlaw_C11(i)- & + constitutive_phenopowerlaw_C12(i)) end select constitutive_phenopowerlaw_Cslip_66(:,:,i) = & math_Mandel3333to66(math_Voigt66to3333(constitutive_phenopowerlaw_Cslip_66(:,:,i))) diff --git a/code/lattice.f90 b/code/lattice.f90 index 078e8c14e..199df9e46 100644 --- a/code/lattice.f90 +++ b/code/lattice.f90 @@ -29,7 +29,7 @@ integer(pInt), parameter :: lattice_maxNinteraction = 20 ! max # of int integer(pInt), pointer, dimension(:,:) :: interactionSlipSlip, & interactionSlipTwin, & interactionTwinSlip, & - interactionTwinTwin + interactionTwinTwin ! Schmid matrices, normal, shear direction and d x n of slip systems real(pReal), allocatable, dimension(:,:,:,:) :: lattice_Sslip @@ -666,10 +666,14 @@ subroutine lattice_init() if(.not. IO_open_file(fileunit,material_configFile)) call IO_error (100) ! corrupt config file Nsections = IO_countSections(fileunit,material_partPhase) -! lattice_Nstructure = 2_pInt + sum(IO_countTagInPart(fileunit,material_partPhase,'covera_ratio',Nsections)) ! fcc + bcc + all hex - lattice_Nstructure = Nsections ! most conservative assumption + lattice_Nstructure = 2_pInt + sum(IO_countTagInPart(fileunit,material_partPhase,'covera_ratio',Nsections)) ! fcc + bcc + all hex +! lattice_Nstructure = Nsections + 2_pInt ! most conservative assumption close(fileunit) + write(6,'(a16,x,i5)') '# sections:',Nsections + write(6,'(a16,x,i5)') '# structures:',lattice_Nstructure + write(6,*) + allocate(lattice_Sslip(3,3,lattice_maxNslip,lattice_Nstructure)); lattice_Sslip = 0.0_pReal allocate(lattice_Sslip_v(6,lattice_maxNslip,lattice_Nstructure)); lattice_Sslip_v = 0.0_pReal allocate(lattice_sd(3,lattice_maxNslip,lattice_Nstructure)); lattice_sd = 0.0_pReal @@ -703,6 +707,7 @@ function lattice_initializeStructure(struct,CoverA) !************************************** use prec, only: pReal,pInt use math + use IO, only: IO_error implicit none character(len=*) struct @@ -818,17 +823,19 @@ function lattice_initializeStructure(struct,CoverA) interactionTwinSlip => lattice_hex_interactionTwinSlip interactionTwinTwin => lattice_hex_interactionTwinTwin endif - end select + end select if (processMe) then - do i = 1,myNslip ! store slip system vectors and Schmid matrix for my structure + if (myStructure > lattice_Nstructure) & + call IO_error(666,0,0,0,'structure index too large') ! check for memory leakage + do i = 1,myNslip ! store slip system vectors and Schmid matrix for my structure lattice_sd(:,i,myStructure) = sd(:,i) lattice_st(:,i,myStructure) = st(:,i) lattice_sn(:,i,myStructure) = sn(:,i) lattice_Sslip(:,:,i,myStructure) = math_tensorproduct(sd(:,i),sn(:,i)) lattice_Sslip_v(:,i,myStructure) = math_Mandel33to6(math_symmetric3x3(lattice_Sslip(:,:,i,myStructure))) enddo - do i = 1,myNtwin ! store twin system vectors and Schmid plus rotation matrix for my structure + do i = 1,myNtwin ! store twin system vectors and Schmid plus rotation matrix for my structure lattice_td(:,i,myStructure) = td(:,i) lattice_tt(:,i,myStructure) = tt(:,i) lattice_tn(:,i,myStructure) = tn(:,i)