diff --git a/code/DAMASK_marc.f90 b/code/DAMASK_marc.f90 index bdeb60964..5a76d9df1 100644 --- a/code/DAMASK_marc.f90 +++ b/code/DAMASK_marc.f90 @@ -1,4 +1,4 @@ -! Copyright 2011 Max-Planck-Institut f√ºr Eisenforschung GmbH +! Copyright 2011 Max-Planck-Institut für Eisenforschung GmbH ! ! This file is part of DAMASK, ! the D√ºsseldorf Advanced MAterial Simulation Kit. diff --git a/code/DAMASK_spectral.f90 b/code/DAMASK_spectral.f90 index c4716b62f..69e71866d 100644 --- a/code/DAMASK_spectral.f90 +++ b/code/DAMASK_spectral.f90 @@ -132,7 +132,7 @@ program DAMASK_spectral 0.0_pReal, 0.0_pReal, 0.0_pReal, 0.0_pReal, 0.0_pReal, 1.0_pReal],& [ 6_pInt, 6_pInt]) real(pReal), dimension(3,3,3,3) :: temp_3333 = 0.0_pReal - integer(pInt) :: size_reduced = 0.0_pReal ! number of stress BCs + integer(pInt) :: size_reduced = 0_pInt ! number of stress BCs !-------------------------------------------------------------------------------------------------- ! pointwise data diff --git a/code/DAMASK_spectral_interface.f90 b/code/DAMASK_spectral_interface.f90 index 219830c6e..deda85309 100644 --- a/code/DAMASK_spectral_interface.f90 +++ b/code/DAMASK_spectral_interface.f90 @@ -35,7 +35,7 @@ CONTAINS ! !******************************************************************** subroutine DAMASK_interface_init() - + use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) implicit none character(len=1024) commandLine, hostName, userName diff --git a/code/constitutive.f90 b/code/constitutive.f90 index 58d6946aa..4db8d12d0 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -101,7 +101,7 @@ logical knownConstitution ! --- PARSE CONSTITUTIONS FROM CONFIG FILE --- if (.not. IO_open_jobFile(fileunit,material_localFileExt)) then ! no local material configuration present... - if (.not. IO_open_file(fileunit,material_configFile)) call IO_error(100) ! ...and cannot open material.config file + if (.not. IO_open_file(fileunit,material_configFile)) call IO_error(100_pInt) ! ...and cannot open material.config file endif call constitutive_j2_init(fileunit) call constitutive_phenopowerlaw_init(fileunit) @@ -114,7 +114,7 @@ close(fileunit) ! --- WRITE DESCRIPTION FILE FOR CONSTITUTIVE PHASE OUTPUT --- if(.not. IO_write_jobFile(fileunit,'outputConstitutive')) then ! problems in writing file - call IO_error (50) + call IO_error (50_pInt) endif do p = 1,material_Nphase i = phase_constitutionInstance(p) ! which instance of a constitution is present phase @@ -323,7 +323,7 @@ endif constitutive_sizePostResults(g,i,e) = constitutive_nonlocal_sizePostResults(myInstance) case default - call IO_error(200,material_phase(g,i,e)) ! unknown constitution + call IO_error(200_pInt,material_phase(g,i,e)) ! unknown constitution end select constitutive_partionedState0(g,i,e)%p = constitutive_state0(g,i,e)%p diff --git a/code/constitutive_dislotwin.f90 b/code/constitutive_dislotwin.f90 index 0d85941e7..f0b651aa8 100644 --- a/code/constitutive_dislotwin.f90 +++ b/code/constitutive_dislotwin.f90 @@ -287,148 +287,148 @@ do ! read thru sections of if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif - if (section > 0 .and. phase_constitution(section) == constitutive_dislotwin_label) then ! one of my sections + if (section > 0_pInt .and. phase_constitution(section) == constitutive_dislotwin_label) then ! one of my sections i = phase_constitutionInstance(section) ! which instance of my constitution is present phase positions = IO_stringPos(line,maxNchunks) tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key select case(tag) case ('(output)') - output = output + 1 - constitutive_dislotwin_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + constitutive_dislotwin_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('lattice_structure') - constitutive_dislotwin_structureName(i) = IO_lc(IO_stringValue(line,positions,2)) + constitutive_dislotwin_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('covera_ratio') - constitutive_dislotwin_CoverA(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_CoverA(i) = IO_floatValue(line,positions,2_pInt) case ('c11') - constitutive_dislotwin_C11(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_C11(i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_dislotwin_C12(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_C12(i) = IO_floatValue(line,positions,2_pInt) case ('c13') - constitutive_dislotwin_C13(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_C13(i) = IO_floatValue(line,positions,2_pInt) case ('c33') - constitutive_dislotwin_C33(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_C33(i) = IO_floatValue(line,positions,2_pInt) case ('c44') - constitutive_dislotwin_C44(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_C44(i) = IO_floatValue(line,positions,2_pInt) case ('nslip') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_Nslip(j,i) = IO_intValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_Nslip(j,i) = IO_intValue(line,positions,1_pInt+j) case ('ntwin') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_dislotwin_Ntwin(j,i) = IO_intValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_dislotwin_Ntwin(j,i) = IO_intValue(line,positions,1_pInt+j) case ('rhoedge0') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_rhoEdge0(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_rhoEdge0(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('rhoedgedip0') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_rhoEdgeDip0(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_rhoEdgeDip0(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('slipburgers') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_burgersPerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_burgersPerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('twinburgers') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_dislotwin_burgersPerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_dislotwin_burgersPerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('qedge') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_QedgePerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_QedgePerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('v0') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_v0PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_v0PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('ndot0') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_dislotwin_Ndot0PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_dislotwin_Ndot0PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('twinsize') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_dislotwin_twinsizePerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_dislotwin_twinsizePerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('clambdaslip') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_dislotwin_CLambdaSlipPerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_dislotwin_CLambdaSlipPerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('grainsize') - constitutive_dislotwin_GrainSize(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_GrainSize(i) = IO_floatValue(line,positions,2_pInt) case ('maxtwinfraction') - constitutive_dislotwin_MaxTwinFraction(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_MaxTwinFraction(i) = IO_floatValue(line,positions,2_pInt) case ('pexponent') - constitutive_dislotwin_p(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_p(i) = IO_floatValue(line,positions,2_pInt) case ('qexponent') - constitutive_dislotwin_q(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_q(i) = IO_floatValue(line,positions,2_pInt) case ('rexponent') - constitutive_dislotwin_r(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_r(i) = IO_floatValue(line,positions,2_pInt) case ('d0') - constitutive_dislotwin_D0(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_D0(i) = IO_floatValue(line,positions,2_pInt) case ('qsd') - constitutive_dislotwin_Qsd(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_Qsd(i) = IO_floatValue(line,positions,2_pInt) case ('atol_rho') - constitutive_dislotwin_aTolRho(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_aTolRho(i) = IO_floatValue(line,positions,2_pInt) case ('cmfptwin') - constitutive_dislotwin_Cmfptwin(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_Cmfptwin(i) = IO_floatValue(line,positions,2_pInt) case ('cthresholdtwin') - constitutive_dislotwin_Cthresholdtwin(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_Cthresholdtwin(i) = IO_floatValue(line,positions,2_pInt) case ('solidsolutionstrength') - constitutive_dislotwin_SolidSolutionStrength(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_SolidSolutionStrength(i) = IO_floatValue(line,positions,2_pInt) case ('l0') - constitutive_dislotwin_L0(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_L0(i) = IO_floatValue(line,positions,2_pInt) case ('cedgedipmindistance') - constitutive_dislotwin_CEdgeDipMinDistance(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_CEdgeDipMinDistance(i) = IO_floatValue(line,positions,2_pInt) case ('catomicvolume') - constitutive_dislotwin_CAtomicVolume(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_CAtomicVolume(i) = IO_floatValue(line,positions,2_pInt) case ('interactionslipslip') - forall (j = 1:lattice_maxNinteraction) & - constitutive_dislotwin_interactionSlipSlip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_dislotwin_interactionSlipSlip(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('interactionsliptwin') - forall (j = 1:lattice_maxNinteraction) & - constitutive_dislotwin_interactionSlipTwin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_dislotwin_interactionSlipTwin(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('interactiontwinslip') - forall (j = 1:lattice_maxNinteraction) & - constitutive_dislotwin_interactionTwinSlip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_dislotwin_interactionTwinSlip(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('interactiontwintwin') - forall (j = 1:lattice_maxNinteraction) & - constitutive_dislotwin_interactionTwinTwin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_dislotwin_interactionTwinTwin(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('sfe_0k') - constitutive_dislotwin_SFE_0K(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_SFE_0K(i) = IO_floatValue(line,positions,2_pInt) case ('dsfe_dt') - constitutive_dislotwin_dSFE_dT(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_dSFE_dT(i) = IO_floatValue(line,positions,2_pInt) case ('shearbandresistance') - constitutive_dislotwin_sbResistance(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_sbResistance(i) = IO_floatValue(line,positions,2_pInt) case ('shearbandvelocity') - constitutive_dislotwin_sbVelocity(i) = IO_floatValue(line,positions,2) + constitutive_dislotwin_sbVelocity(i) = IO_floatValue(line,positions,2_pInt) end select endif enddo -100 do i = 1,maxNinstance +100 do i = 1_pInt,maxNinstance constitutive_dislotwin_structure(i) = & lattice_initializeStructure(constitutive_dislotwin_structureName(i),constitutive_dislotwin_CoverA(i)) myStructure = constitutive_dislotwin_structure(i) !* Sanity checks - if (myStructure < 1 .or. myStructure > 3) call IO_error(205) - if (sum(constitutive_dislotwin_Nslip(:,i)) <= 0_pInt) call IO_error(225) - if (sum(constitutive_dislotwin_Ntwin(:,i)) < 0_pInt) call IO_error(225) !*** + if (myStructure < 1_pInt .or. myStructure > 3_pInt) call IO_error(205_pInt) + if (sum(constitutive_dislotwin_Nslip(:,i)) <= 0_pInt) call IO_error(225_pInt) + if (sum(constitutive_dislotwin_Ntwin(:,i)) < 0_pInt) call IO_error(225_pInt) !*** do f = 1,lattice_maxNslipFamily if (constitutive_dislotwin_Nslip(f,i) > 0_pInt) then - if (constitutive_dislotwin_rhoEdge0(f,i) < 0.0_pReal) call IO_error(220) - if (constitutive_dislotwin_rhoEdgeDip0(f,i) < 0.0_pReal) call IO_error(220) - if (constitutive_dislotwin_burgersPerSlipFamily(f,i) <= 0.0_pReal) call IO_error(221) - if (constitutive_dislotwin_v0PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(226) + if (constitutive_dislotwin_rhoEdge0(f,i) < 0.0_pReal) call IO_error(220_pInt) + if (constitutive_dislotwin_rhoEdgeDip0(f,i) < 0.0_pReal) call IO_error(220_pInt) + if (constitutive_dislotwin_burgersPerSlipFamily(f,i) <= 0.0_pReal) call IO_error(221_pInt) + if (constitutive_dislotwin_v0PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(226_pInt) endif enddo do f = 1,lattice_maxNtwinFamily if (constitutive_dislotwin_Ntwin(f,i) > 0_pInt) then - if (constitutive_dislotwin_burgersPerTwinFamily(f,i) <= 0.0_pReal) call IO_error(221) !*** - if (constitutive_dislotwin_Ndot0PerTwinFamily(f,i) < 0.0_pReal) call IO_error(226) !*** + if (constitutive_dislotwin_burgersPerTwinFamily(f,i) <= 0.0_pReal) call IO_error(221_pInt) !*** + if (constitutive_dislotwin_Ndot0PerTwinFamily(f,i) < 0.0_pReal) call IO_error(226_pInt) !*** endif enddo ! if (any(constitutive_dislotwin_interactionSlipSlip(1:maxval(lattice_interactionSlipSlip(:,:,myStructure)),i) < 1.0_pReal)) call IO_error(229) - if (constitutive_dislotwin_CAtomicVolume(i) <= 0.0_pReal) call IO_error(230) - if (constitutive_dislotwin_D0(i) <= 0.0_pReal) call IO_error(231) - if (constitutive_dislotwin_Qsd(i) <= 0.0_pReal) call IO_error(232) - if (constitutive_dislotwin_aTolRho(i) <= 0.0_pReal) call IO_error(233) - if (constitutive_dislotwin_sbResistance(i) <= 0.0_pReal) call IO_error(234) - if (constitutive_dislotwin_sbVelocity(i) < 0.0_pReal) call IO_error(235) - if (constitutive_dislotwin_SFE_0K(i) == 0.0_pReal .AND. & - constitutive_dislotwin_dSFE_dT(i) == 0.0_pReal) call IO_error(223) + if (constitutive_dislotwin_CAtomicVolume(i) <= 0.0_pReal) call IO_error(230_pInt) + if (constitutive_dislotwin_D0(i) <= 0.0_pReal) call IO_error(231_pInt) + if (constitutive_dislotwin_Qsd(i) <= 0.0_pReal) call IO_error(232_pInt) + if (constitutive_dislotwin_aTolRho(i) <= 0.0_pReal) call IO_error(233_pInt) + if (constitutive_dislotwin_sbResistance(i) <= 0.0_pReal) call IO_error(234_pInt) + if (constitutive_dislotwin_sbVelocity(i) < 0.0_pReal) call IO_error(235_pInt) + if (constitutive_dislotwin_SFE_0K(i) == 0.0_pReal .and. & + constitutive_dislotwin_dSFE_dT(i) == 0.0_pReal) call IO_error(223_pInt) !* Determine total number of active slip or twin systems constitutive_dislotwin_Nslip(:,i) = min(lattice_NslipSystem(:,myStructure),constitutive_dislotwin_Nslip(:,i)) @@ -473,25 +473,25 @@ allocate(constitutive_dislotwin_Ctwin_3333(3,3,3,3,maxTotalNtwin,maxNinstance)) constitutive_dislotwin_Ctwin_66 = 0.0_pReal constitutive_dislotwin_Ctwin_3333 = 0.0_pReal -do i = 1,maxNinstance +do i = 1_pInt,maxNinstance myStructure = constitutive_dislotwin_structure(i) !* Inverse lookup of my slip system family l = 0_pInt - do f = 1,lattice_maxNslipFamily - do k = 1,constitutive_dislotwin_Nslip(f,i) - l = l + 1 + do f = 1_pInt,lattice_maxNslipFamily + do k = 1_pInt,constitutive_dislotwin_Nslip(f,i) + l = l + 1_pInt constitutive_dislotwin_slipFamily(l,i) = f - constitutive_dislotwin_slipSystemLattice(l,i) = sum(lattice_NslipSystem(1:f-1,myStructure)) + k + constitutive_dislotwin_slipSystemLattice(l,i) = sum(lattice_NslipSystem(1:f-1_pInt,myStructure)) + k enddo; enddo !* Inverse lookup of my twin system family l = 0_pInt - do f = 1,lattice_maxNtwinFamily - do k = 1,constitutive_dislotwin_Ntwin(f,i) - l = l + 1 + do f = 1_pInt,lattice_maxNtwinFamily + do k = 1_pInt,constitutive_dislotwin_Ntwin(f,i) + l = l + 1_pInt constitutive_dislotwin_twinFamily(l,i) = f - constitutive_dislotwin_twinSystemLattice(l,i) = sum(lattice_NtwinSystem(1:f-1,myStructure)) + k + constitutive_dislotwin_twinSystemLattice(l,i) = sum(lattice_NtwinSystem(1:f-1_pInt,myStructure)) + k enddo; enddo !* Determine size of state array @@ -504,7 +504,7 @@ do i = 1,maxNinstance size(constitutive_dislotwin_listDependentSlipStates)*ns+size(constitutive_dislotwin_listDependentTwinStates)*nt !* Determine size of postResults array - do o = 1,maxval(phase_Noutput) + do o = 1_pInt,maxval(phase_Noutput) select case(constitutive_dislotwin_output(o,i)) case('edge_density', & 'dipole_density', & @@ -545,14 +545,14 @@ do i = 1,maxNinstance !* Elasticity matrix and shear modulus according to material.config select case (myStructure) - case(1:2) ! cubic(s) - forall(k=1:3) - forall(j=1:3) & + case(1_pInt:2_pInt) ! cubic(s) + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) & constitutive_dislotwin_Cslip_66(k,j,i) = constitutive_dislotwin_C12(i) constitutive_dislotwin_Cslip_66(k,k,i) = constitutive_dislotwin_C11(i) - constitutive_dislotwin_Cslip_66(k+3,k+3,i) = constitutive_dislotwin_C44(i) + constitutive_dislotwin_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_dislotwin_C44(i) end forall - case(3:) ! all hex + case(3_pInt:) ! all hex constitutive_dislotwin_Cslip_66(1,1,i) = constitutive_dislotwin_C11(i) constitutive_dislotwin_Cslip_66(2,2,i) = constitutive_dislotwin_C11(i) constitutive_dislotwin_Cslip_66(3,3,i) = constitutive_dislotwin_C33(i) @@ -572,23 +572,24 @@ do i = 1,maxNinstance 0.2_pReal*(constitutive_dislotwin_C11(i)-constitutive_dislotwin_C12(i))+0.3_pReal*constitutive_dislotwin_C44(i) !* Construction of the twin elasticity matrices - do j=1,lattice_maxNtwinFamily - do k=1,constitutive_dislotwin_Ntwin(j,i) - do l=1,3 ; do m=1,3 ; do n=1,3 ; do o=1,3 ; do p=1,3 ; do q=1,3 ; do r=1,3 ; do s=1,3 - constitutive_dislotwin_Ctwin_3333(l,m,n,o,sum(constitutive_dislotwin_Nslip(1:j-1,i))+k,i) = & - constitutive_dislotwin_Ctwin_3333(l,m,n,o,sum(constitutive_dislotwin_Nslip(1:j-1,i))+k,i) + & - constitutive_dislotwin_Cslip_3333(p,q,r,s,i)*& - lattice_Qtwin(l,p,sum(lattice_NslipSystem(1:j-1,myStructure))+k,myStructure)* & - lattice_Qtwin(m,q,sum(lattice_NslipSystem(1:j-1,myStructure))+k,myStructure)* & - lattice_Qtwin(n,r,sum(lattice_NslipSystem(1:j-1,myStructure))+k,myStructure)* & - lattice_Qtwin(o,s,sum(lattice_NslipSystem(1:j-1,myStructure))+k,myStructure) + do j=1_pInt,lattice_maxNtwinFamily + do k=1_pInt,constitutive_dislotwin_Ntwin(j,i) + do l=1_pInt,3_pInt ; do m=1_pInt,3_pInt ; do n=1_pInt,3_pInt ; do o=1_pInt,3_pInt + do p=1_pInt,3_pInt ; do q=1_pInt,3_pInt ; do r=1_pInt,3_pInt ; do s=1_pInt,3_pInt + constitutive_dislotwin_Ctwin_3333(l,m,n,o,sum(constitutive_dislotwin_Nslip(1:j-1_pInt,i))+k,i) = & + constitutive_dislotwin_Ctwin_3333(l,m,n,o,sum(constitutive_dislotwin_Nslip(1:j-1_pInt,i))+k,i) + & + constitutive_dislotwin_Cslip_3333(p,q,r,s,i)*& + lattice_Qtwin(l,p,sum(lattice_NslipSystem(1:j-1_pInt,myStructure))+k,myStructure)* & + lattice_Qtwin(m,q,sum(lattice_NslipSystem(1:j-1_pInt,myStructure))+k,myStructure)* & + lattice_Qtwin(n,r,sum(lattice_NslipSystem(1:j-1_pInt,myStructure))+k,myStructure)* & + lattice_Qtwin(o,s,sum(lattice_NslipSystem(1:j-1_pInt,myStructure))+k,myStructure) enddo ; enddo ; enddo ; enddo ; enddo ; enddo ; enddo ; enddo constitutive_dislotwin_Ctwin_66(:,:,k,i) = math_Mandel3333to66(constitutive_dislotwin_Ctwin_3333(:,:,:,:,k,i)) enddo enddo !* Burgers vector, dislocation velocity prefactor, mean free path prefactor and minimum dipole distance for each slip system - do s = 1,constitutive_dislotwin_totalNslip(i) + do s = 1_pInt,constitutive_dislotwin_totalNslip(i) f = constitutive_dislotwin_slipFamily(s,i) constitutive_dislotwin_burgersPerSlipSystem(s,i) = constitutive_dislotwin_burgersPerSlipFamily(f,i) constitutive_dislotwin_QedgePerSlipSystem(s,i) = constitutive_dislotwin_QedgePerSlipFamily(f,i) @@ -597,7 +598,7 @@ do i = 1,maxNinstance enddo !* Burgers vector, nucleation rate prefactor and twin size for each twin system - do s = 1,constitutive_dislotwin_totalNtwin(i) + do s = 1_pInt,constitutive_dislotwin_totalNtwin(i) f = constitutive_dislotwin_twinFamily(s,i) constitutive_dislotwin_burgersPerTwinSystem(s,i) = constitutive_dislotwin_burgersPerTwinFamily(f,i) constitutive_dislotwin_Ndot0PerTwinSystem(s,i) = constitutive_dislotwin_Ndot0PerTwinFamily(f,i) @@ -605,41 +606,41 @@ do i = 1,maxNinstance enddo !* Construction of interaction matrices - do s1 = 1,constitutive_dislotwin_totalNslip(i) - do s2 = 1,constitutive_dislotwin_totalNslip(i) + do s1 = 1_pInt,constitutive_dislotwin_totalNslip(i) + do s2 = 1_pInt,constitutive_dislotwin_totalNslip(i) constitutive_dislotwin_interactionMatrixSlipSlip(s1,s2,i) = & constitutive_dislotwin_interactionSlipSlip(lattice_interactionSlipSlip(constitutive_dislotwin_slipSystemLattice(s1,i), & constitutive_dislotwin_slipSystemLattice(s2,i), & myStructure),i) enddo; enddo - do s1 = 1,constitutive_dislotwin_totalNslip(i) - do t2 = 1,constitutive_dislotwin_totalNtwin(i) + do s1 = 1_pInt,constitutive_dislotwin_totalNslip(i) + do t2 = 1_pInt,constitutive_dislotwin_totalNtwin(i) constitutive_dislotwin_interactionMatrixSlipTwin(s1,t2,i) = & - constitutive_dislotwin_interactionSlipTwin(lattice_interactionSlipTwin(constitutive_dislotwin_slipSystemLattice(s1,i), & - constitutive_dislotwin_twinSystemLattice(t2,i), & - myStructure),i) + constitutive_dislotwin_interactionSlipTwin(& + lattice_interactionSlipTwin(constitutive_dislotwin_slipSystemLattice(s1,i), & + constitutive_dislotwin_twinSystemLattice(t2,i),myStructure),i) enddo; enddo - do t1 = 1,constitutive_dislotwin_totalNtwin(i) - do s2 = 1,constitutive_dislotwin_totalNslip(i) + do t1 = 1_pInt,constitutive_dislotwin_totalNtwin(i) + do s2 = 1_pInt,constitutive_dislotwin_totalNslip(i) constitutive_dislotwin_interactionMatrixTwinSlip(t1,s2,i) = & - constitutive_dislotwin_interactionTwinSlip(lattice_interactionTwinSlip(constitutive_dislotwin_twinSystemLattice(t1,i), & - constitutive_dislotwin_slipSystemLattice(s2,i), & - myStructure),i) + constitutive_dislotwin_interactionTwinSlip(lattice_interactionTwinSlip(& + constitutive_dislotwin_twinSystemLattice(t1,i), & + constitutive_dislotwin_slipSystemLattice(s2,i), myStructure),i) enddo; enddo - do t1 = 1,constitutive_dislotwin_totalNtwin(i) - do t2 = 1,constitutive_dislotwin_totalNtwin(i) + do t1 = 1_pInt,constitutive_dislotwin_totalNtwin(i) + do t2 = 1_pInt,constitutive_dislotwin_totalNtwin(i) constitutive_dislotwin_interactionMatrixTwinTwin(t1,t2,i) = & - constitutive_dislotwin_interactionTwinTwin(lattice_interactionTwinTwin(constitutive_dislotwin_twinSystemLattice(t1,i), & - constitutive_dislotwin_twinSystemLattice(t2,i), & - myStructure),i) + constitutive_dislotwin_interactionTwinTwin(& + lattice_interactionTwinTwin(constitutive_dislotwin_twinSystemLattice(t1,i), & + constitutive_dislotwin_twinSystemLattice(t2,i), myStructure),i) enddo; enddo !* Calculation of forest projections for edge dislocations - do s1 = 1,constitutive_dislotwin_totalNslip(i) - do s2 = 1,constitutive_dislotwin_totalNslip(i) + do s1 = 1_pInt,constitutive_dislotwin_totalNslip(i) + do s2 = 1_pInt,constitutive_dislotwin_totalNslip(i) constitutive_dislotwin_forestProjectionEdge(s1,s2,i) = & abs(math_mul3x3(lattice_sn(:,constitutive_dislotwin_slipSystemLattice(s1,i),myStructure), & lattice_st(:,constitutive_dislotwin_slipSystemLattice(s2,i),myStructure))) @@ -678,7 +679,7 @@ constitutive_dislotwin_stateInit = 0.0_pReal !* Initialize basic slip state variables s1 = 0_pInt -do f = 1,lattice_maxNslipFamily +do f = 1_pInt,lattice_maxNslipFamily s0 = s1 + 1_pInt s1 = s0 + constitutive_dislotwin_Nslip(f,myInstance) - 1_pInt do s = s0,s1 @@ -686,25 +687,25 @@ do f = 1,lattice_maxNslipFamily rhoEdgeDip0(s) = constitutive_dislotwin_rhoEdgeDip0(f,myInstance) enddo enddo -constitutive_dislotwin_stateInit(1:ns) = rhoEdge0 -constitutive_dislotwin_stateInit(ns+1:2*ns) = rhoEdgeDip0 +constitutive_dislotwin_stateInit(1:ns) = rhoEdge0 +constitutive_dislotwin_stateInit(ns+1:2_pInt*ns) = rhoEdgeDip0 !* Initialize dependent slip microstructural variables -forall (s = 1:ns) & +forall (s = 1_pInt:ns) & invLambdaSlip0(s) = sqrt(dot_product((rhoEdge0+rhoEdgeDip0),constitutive_dislotwin_forestProjectionEdge(1:ns,s,myInstance)))/ & constitutive_dislotwin_CLambdaSlipPerSlipSystem(s,myInstance) -constitutive_dislotwin_stateInit(2*ns+nt+1:3*ns+nt) = invLambdaSlip0 +constitutive_dislotwin_stateInit(2_pInt*ns+nt+1:3_pInt*ns+nt) = invLambdaSlip0 -forall (s = 1:ns) & +forall (s = 1_pInt:ns) & MeanFreePathSlip0(s) = & constitutive_dislotwin_GrainSize(myInstance)/(1.0_pReal+invLambdaSlip0(s)*constitutive_dislotwin_GrainSize(myInstance)) -constitutive_dislotwin_stateInit(4*ns+2*nt+1:5*ns+2*nt) = MeanFreePathSlip0 +constitutive_dislotwin_stateInit(4_pInt*ns+2_pInt*nt+1:5_pInt*ns+2_pInt*nt) = MeanFreePathSlip0 -forall (s = 1:ns) & +forall (s = 1_pInt:ns) & tauSlipThreshold0(s) = constitutive_dislotwin_SolidSolutionStrength(myInstance)+ & constitutive_dislotwin_Gmod(myInstance)*constitutive_dislotwin_burgersPerSlipSystem(s,myInstance)* & sqrt(dot_product((rhoEdge0+rhoEdgeDip0),constitutive_dislotwin_interactionMatrixSlipSlip(1:ns,s,myInstance))) -constitutive_dislotwin_stateInit(5*ns+3*nt+1:6*ns+3*nt) = tauSlipThreshold0 +constitutive_dislotwin_stateInit(5_pInt*ns+3_pInt*nt+1:6_pInt*ns+3_pInt*nt) = tauSlipThreshold0 !* Initialize dependent twin microstructural variables forall (t = 1:nt) & diff --git a/code/constitutive_j2.f90 b/code/constitutive_j2.f90 index 8093f6c0d..127451dab 100644 --- a/code/constitutive_j2.f90 +++ b/code/constitutive_j2.f90 @@ -83,14 +83,15 @@ subroutine constitutive_j2_init(file) !************************************** !* Module initialization * !************************************** + use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) use prec, only: pInt, pReal use math, only: math_Mandel3333to66, math_Voigt66to3333 use IO use material use debug, only: debug_verbosity integer(pInt), intent(in) :: file - integer(pInt), parameter :: maxNchunks = 7 - integer(pInt), dimension(1+2*maxNchunks) :: positions + integer(pInt), parameter :: maxNchunks = 7_pInt + integer(pInt), dimension(1+2_pInt*maxNchunks) :: positions integer(pInt) section, maxNinstance, i,j,k, output, mySize character(len=64) tag character(len=1024) line @@ -103,9 +104,9 @@ subroutine constitutive_j2_init(file) !$OMP END CRITICAL (write2out) maxNinstance = count(phase_constitution == constitutive_j2_label) - if (maxNinstance == 0) return + if (maxNinstance == 0_pInt) return - if (debug_verbosity > 0) then + if (debug_verbosity > 0_pInt) then !$OMP CRITICAL (write2out) write(6,'(a16,1x,i5)') '# instances:',maxNinstance write(6,*) @@ -126,12 +127,12 @@ subroutine constitutive_j2_init(file) allocate(constitutive_j2_n(maxNinstance)) ; constitutive_j2_n = 0.0_pReal allocate(constitutive_j2_h0(maxNinstance)) ; constitutive_j2_h0 = 0.0_pReal allocate(constitutive_j2_tausat(maxNinstance)) ; constitutive_j2_tausat = 0.0_pReal - allocate(constitutive_j2_a(maxNinstance)) ; constitutive_j2_a = 0.0_pReal + allocate(constitutive_j2_a(maxNinstance)) ; constitutive_j2_a = 0.0_pReal allocate(constitutive_j2_aTolResistance(maxNinstance)) ; constitutive_j2_aTolResistance = 0.0_pReal rewind(file) line = '' - section = 0 + section = 0_pInt do while (IO_lc(IO_getTag(line,'<','>')) /= 'phase') ! wind forward to read(file,'(a1024)',END=100) line @@ -142,53 +143,53 @@ subroutine constitutive_j2_init(file) if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif - if (section > 0 .and. phase_constitution(section) == constitutive_j2_label) then ! one of my sections + if (section > 0_pInt .and. phase_constitution(section) == constitutive_j2_label) then ! one of my sections i = phase_constitutionInstance(section) ! which instance of my constitution is present phase positions = IO_stringPos(line,maxNchunks) tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key select case(tag) case ('(output)') - output = output + 1 - constitutive_j2_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + constitutive_j2_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('c11') - constitutive_j2_C11(i) = IO_floatValue(line,positions,2) + constitutive_j2_C11(i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_j2_C12(i) = IO_floatValue(line,positions,2) + constitutive_j2_C12(i) = IO_floatValue(line,positions,2_pInt) case ('tau0') - constitutive_j2_tau0(i) = IO_floatValue(line,positions,2) + constitutive_j2_tau0(i) = IO_floatValue(line,positions,2_pInt) case ('gdot0') - constitutive_j2_gdot0(i) = IO_floatValue(line,positions,2) + constitutive_j2_gdot0(i) = IO_floatValue(line,positions,2_pInt) case ('n') - constitutive_j2_n(i) = IO_floatValue(line,positions,2) + constitutive_j2_n(i) = IO_floatValue(line,positions,2_pInt) case ('h0') - constitutive_j2_h0(i) = IO_floatValue(line,positions,2) + constitutive_j2_h0(i) = IO_floatValue(line,positions,2_pInt) case ('tausat') - constitutive_j2_tausat(i) = IO_floatValue(line,positions,2) + constitutive_j2_tausat(i) = IO_floatValue(line,positions,2_pInt) case ('a', 'w0') - constitutive_j2_a(i) = IO_floatValue(line,positions,2) + constitutive_j2_a(i) = IO_floatValue(line,positions,2_pInt) case ('taylorfactor') - constitutive_j2_fTaylor(i) = IO_floatValue(line,positions,2) + constitutive_j2_fTaylor(i) = IO_floatValue(line,positions,2_pInt) case ('atol_resistance') - constitutive_j2_aTolResistance(i) = IO_floatValue(line,positions,2) + constitutive_j2_aTolResistance(i) = IO_floatValue(line,positions,2_pInt) end select endif enddo -100 do i = 1,maxNinstance ! sanity checks - if (constitutive_j2_tau0(i) < 0.0_pReal) call IO_error(210) - if (constitutive_j2_gdot0(i) <= 0.0_pReal) call IO_error(211) - if (constitutive_j2_n(i) <= 0.0_pReal) call IO_error(212) - if (constitutive_j2_tausat(i) <= 0.0_pReal) call IO_error(213) - if (constitutive_j2_a(i) <= 0.0_pReal) call IO_error(241) - if (constitutive_j2_fTaylor(i) <= 0.0_pReal) call IO_error(240) - if (constitutive_j2_aTolResistance(i) <= 0.0_pReal) call IO_error(242) +100 do i = 1_pInt,maxNinstance ! sanity checks + if (constitutive_j2_tau0(i) < 0.0_pReal) call IO_error(210_pInt) + if (constitutive_j2_gdot0(i) <= 0.0_pReal) call IO_error(211_pInt) + if (constitutive_j2_n(i) <= 0.0_pReal) call IO_error(212_pInt) + if (constitutive_j2_tausat(i) <= 0.0_pReal) call IO_error(213_pInt) + if (constitutive_j2_a(i) <= 0.0_pReal) call IO_error(241_pInt) + if (constitutive_j2_fTaylor(i) <= 0.0_pReal) call IO_error(240_pInt) + if (constitutive_j2_aTolResistance(i) <= 0.0_pReal) call IO_error(242_pInt) enddo - do i = 1,maxNinstance - do j = 1,maxval(phase_Noutput) + do i = 1_pInt,maxNinstance + do j = 1_pInt,maxval(phase_Noutput) select case(constitutive_j2_output(j,i)) case('flowstress') mySize = 1_pInt @@ -205,22 +206,24 @@ subroutine constitutive_j2_init(file) endif enddo - constitutive_j2_sizeDotState(i) = 1 - constitutive_j2_sizeState(i) = 1 + constitutive_j2_sizeDotState(i) = 1_pInt + constitutive_j2_sizeState(i) = 1_pInt - forall(k=1:3) - forall(j=1:3) & + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) & constitutive_j2_Cslip_66(k,j,i) = constitutive_j2_C12(i) constitutive_j2_Cslip_66(k,k,i) = constitutive_j2_C11(i) - constitutive_j2_Cslip_66(k+3,k+3,i) = 0.5_pReal*(constitutive_j2_C11(i)-constitutive_j2_C12(i)) + constitutive_j2_Cslip_66(k+3_pInt,k+3_pInt,i) = 0.5_pReal*(constitutive_j2_C11(i)-constitutive_j2_C12(i)) end forall constitutive_j2_Cslip_66(1:6,1:6,i) = & math_Mandel3333to66(math_Voigt66to3333(constitutive_j2_Cslip_66(1:6,1:6,i))) + enddo return + endsubroutine @@ -241,6 +244,7 @@ pure function constitutive_j2_stateInit(myInstance) endfunction + !********************************************************************* !* relevant microstructural state * !********************************************************************* @@ -261,6 +265,7 @@ real(pReal), dimension(constitutive_j2_sizeState(myInstance)) :: & constitutive_j2_aTolState = constitutive_j2_aTolResistance(myInstance) + endfunction @@ -288,6 +293,7 @@ function constitutive_j2_homogenizedC(state,ipc,ip,el) return + endfunction @@ -310,6 +316,7 @@ subroutine constitutive_j2_microstructure(Temperature,state,ipc,ip,el) real(pReal) Temperature type(p_vec), dimension(homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems) :: state + matID = phase_constitutionInstance(material_phase(ipc,ip,el)) endsubroutine @@ -335,6 +342,7 @@ pure subroutine constitutive_j2_LpAndItsTangent(Lp, dLp_dTstar_99, Tstar_dev_v, material_phase, & phase_constitutionInstance + implicit none !*** input variables ***! @@ -458,6 +466,7 @@ pure function constitutive_j2_dotState(Tstar_v, Temperature, state, g, ip, el) endfunction + !**************************************************************** !* calculates the rate of change of temperature * !**************************************************************** diff --git a/code/constitutive_nonlocal.f90 b/code/constitutive_nonlocal.f90 index c6d4ac69c..c6b87d0dd 100644 --- a/code/constitutive_nonlocal.f90 +++ b/code/constitutive_nonlocal.f90 @@ -140,6 +140,7 @@ CONTAINS !************************************** subroutine constitutive_nonlocal_init(file) +use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) use prec, only: pInt, pReal use math, only: math_Mandel3333to66, & math_Voigt66to3333, & @@ -328,7 +329,7 @@ constitutive_nonlocal_peierlsStressPerSlipFamily = 0.0_pReal rewind(file) line = '' -section = 0 +section = 0_pInt do while (IO_lc(IO_getTag(line,'<','>')) /= 'phase') ! wind forward to read(file,'(a1024)',END=100) line @@ -339,98 +340,108 @@ do if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter cycle endif - if (section > 0 .and. phase_constitution(section) == constitutive_nonlocal_label) then ! one of my sections + if (section > 0_pInt .and. phase_constitution(section) == constitutive_nonlocal_label) then ! one of my sections i = phase_constitutionInstance(section) ! which instance of my constitution is present phase positions = IO_stringPos(line,maxNchunks) - tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key select case(tag) case('constitution','/nonlocal/') cycle case ('(output)') - output = output + 1 - constitutive_nonlocal_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + constitutive_nonlocal_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('lattice_structure') - constitutive_nonlocal_structureName(i) = IO_lc(IO_stringValue(line,positions,2)) + constitutive_nonlocal_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('c/a_ratio','covera_ratio') - constitutive_nonlocal_CoverA(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_CoverA(i) = IO_floatValue(line,positions,2_pInt) case ('c11') - constitutive_nonlocal_C11(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_C11(i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_nonlocal_C12(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_C12(i) = IO_floatValue(line,positions,2_pInt) case ('c13') - constitutive_nonlocal_C13(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_C13(i) = IO_floatValue(line,positions,2_pInt) case ('c33') - constitutive_nonlocal_C33(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_C33(i) = IO_floatValue(line,positions,2_pInt) case ('c44') - constitutive_nonlocal_C44(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_C44(i) = IO_floatValue(line,positions,2_pInt) case ('nslip') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_Nslip(f,i) = IO_intValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_Nslip(f,i)& + = IO_intValue(line,positions,1_pInt+f) case ('rhosgledgepos0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglEdgePos0(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglEdgePos0(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('rhosgledgeneg0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglEdgeNeg0(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglEdgeNeg0(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('rhosglscrewpos0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglScrewPos0(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglScrewPos0(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('rhosglscrewneg0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglScrewNeg0(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_rhoSglScrewNeg0(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('rhodipedge0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_rhoDipEdge0(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_rhoDipEdge0(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('rhodipscrew0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_rhoDipScrew0(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_rhoDipScrew0(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('lambda0') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_lambda0PerSlipFamily(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_lambda0PerSlipFamily(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case ('burgers') - forall (f = 1:lattice_maxNslipFamily) constitutive_nonlocal_burgersPerSlipFamily(f,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_burgersPerSlipFamily(f,i)& + = IO_floatValue(line,positions,1_pInt+f) case('cutoffradius','r') - constitutive_nonlocal_R(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_R(i) = IO_floatValue(line,positions,2_pInt) case('minimumdipoleheightedge','ddipminedge') - forall (f = 1:lattice_maxNslipFamily) & - constitutive_nonlocal_minimumDipoleHeightPerSlipFamily(f,1,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_minimumDipoleHeightPerSlipFamily(f,1_pInt,i)& + = IO_floatValue(line,positions,1_pInt+f) case('minimumdipoleheightscrew','ddipminscrew') - forall (f = 1:lattice_maxNslipFamily) & - constitutive_nonlocal_minimumDipoleHeightPerSlipFamily(f,2,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_minimumDipoleHeightPerSlipFamily(f,2_pInt,i)& + = IO_floatValue(line,positions,1_pInt+f) case('atomicvolume') - constitutive_nonlocal_atomicVolume(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_atomicVolume(i) = IO_floatValue(line,positions,2_pInt) case('selfdiffusionprefactor','dsd0') - constitutive_nonlocal_Dsd0(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_Dsd0(i) = IO_floatValue(line,positions,2_pInt) case('selfdiffusionenergy','qsd') - constitutive_nonlocal_Qsd(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_Qsd(i) = IO_floatValue(line,positions,2_pInt) case('atol_rho') - constitutive_nonlocal_aTolRho(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_aTolRho(i) = IO_floatValue(line,positions,2_pInt) case ('interaction_slipslip') - forall (it = 1:lattice_maxNinteraction) constitutive_nonlocal_interactionSlipSlip(it,i) = IO_floatValue(line,positions,1+it) + forall (it = 1_pInt:lattice_maxNinteraction) constitutive_nonlocal_interactionSlipSlip(it,i)& + = IO_floatValue(line,positions,1_pInt+it) case('peierlsstressedge') - forall (f = 1:lattice_maxNslipFamily) & - constitutive_nonlocal_peierlsStressPerSlipFamily(f,1,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_peierlsStressPerSlipFamily(f,1_pInt,i)& + = IO_floatValue(line,positions,1_pInt+f) case('peierlsstressscrew') - forall (f = 1:lattice_maxNslipFamily) & - constitutive_nonlocal_peierlsStressPerSlipFamily(f,2,i) = IO_floatValue(line,positions,1+f) + forall (f = 1_pInt:lattice_maxNslipFamily) constitutive_nonlocal_peierlsStressPerSlipFamily(f,2_pInt,i)& + = IO_floatValue(line,positions,1_pInt+f) case('doublekinkwidth') - constitutive_nonlocal_doublekinkwidth(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_doublekinkwidth(i) = IO_floatValue(line,positions,2_pInt) case('solidsolutionenergy') - constitutive_nonlocal_solidSolutionEnergy(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_solidSolutionEnergy(i) = IO_floatValue(line,positions,2_pInt) case('solidsolutionsize') - constitutive_nonlocal_solidSolutionSize(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_solidSolutionSize(i) = IO_floatValue(line,positions,2_pInt) case('solidsolutionconcentration') - constitutive_nonlocal_solidSolutionConcentration(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_solidSolutionConcentration(i) = IO_floatValue(line,positions,2_pInt) case('p') - constitutive_nonlocal_p(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_p(i) = IO_floatValue(line,positions,2_pInt) case('q') - constitutive_nonlocal_q(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_q(i) = IO_floatValue(line,positions,2_pInt) case('viscosity','glideviscosity') - constitutive_nonlocal_viscosity(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_viscosity(i) = IO_floatValue(line,positions,2_pInt) case('attackfrequency','fattack') - constitutive_nonlocal_fattack(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_fattack(i) = IO_floatValue(line,positions,2_pInt) case('rhosglscatter') - constitutive_nonlocal_rhoSglScatter(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_rhoSglScatter(i) = IO_floatValue(line,positions,2_pInt) case('surfacetransmissivity') - constitutive_nonlocal_surfaceTransmissivity(i) = IO_floatValue(line,positions,2) + constitutive_nonlocal_surfaceTransmissivity(i) = IO_floatValue(line,positions,2_pInt) case default - call IO_error(236,ext_msg=tag) + call IO_error(236_pInt,ext_msg=tag) end select endif enddo @@ -651,13 +662,13 @@ do i = 1,maxNinstance !*** elasticity matrix and shear modulus according to material.config select case (myStructure) - case(1:2) ! cubic(s) - forall(k=1:3) - forall(j=1:3) constitutive_nonlocal_Cslip_66(k,j,i) = constitutive_nonlocal_C12(i) + case(1_pInt:2_pInt) ! cubic(s) + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) constitutive_nonlocal_Cslip_66(k,j,i) = constitutive_nonlocal_C12(i) constitutive_nonlocal_Cslip_66(k,k,i) = constitutive_nonlocal_C11(i) - constitutive_nonlocal_Cslip_66(k+3,k+3,i) = constitutive_nonlocal_C44(i) + constitutive_nonlocal_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_nonlocal_C44(i) end forall - case(3:) ! all hex + case(3_pInt:) ! all hex constitutive_nonlocal_Cslip_66(1,1,i) = constitutive_nonlocal_C11(i) constitutive_nonlocal_Cslip_66(2,2,i) = constitutive_nonlocal_C11(i) constitutive_nonlocal_Cslip_66(3,3,i) = constitutive_nonlocal_C33(i) @@ -681,7 +692,7 @@ do i = 1,maxNinstance / ( 4.0_pReal*constitutive_nonlocal_C11(i) + 6.0_pReal*constitutive_nonlocal_C12(i) & + 2.0_pReal*constitutive_nonlocal_C44(i) ) ! C12iso/(C11iso+C12iso) with C11iso=(3*C11+2*C12+4*C44)/5 and C12iso=(C11+4*C12-2*C44)/5 - do s1 = 1,ns + do s1 = 1_pInt,ns f = constitutive_nonlocal_slipFamily(s1,i) !*** burgers vector, mean free path prefactor and minimum dipole distance for each slip system @@ -708,17 +719,16 @@ do i = 1,maxNinstance constitutive_nonlocal_interactionMatrixSlipSlip(s1,s2,i) & = constitutive_nonlocal_interactionSlipSlip(lattice_interactionSlipSlip(constitutive_nonlocal_slipSystemLattice(s1,i), & constitutive_nonlocal_slipSystemLattice(s2,i), & - myStructure), & - i) + myStructure), i) enddo !*** rotation matrix from lattice configuration to slip system constitutive_nonlocal_lattice2slip(1:3,1:3,s1,i) & - = math_transpose33( reshape((/ lattice_sd(1:3, constitutive_nonlocal_slipSystemLattice(s1,i), myStructure), & + = math_transpose33( reshape([ lattice_sd(1:3, constitutive_nonlocal_slipSystemLattice(s1,i), myStructure), & -lattice_st(1:3, constitutive_nonlocal_slipSystemLattice(s1,i), myStructure), & - lattice_sn(1:3, constitutive_nonlocal_slipSystemLattice(s1,i), myStructure)/), (/3,3/))) + lattice_sn(1:3, constitutive_nonlocal_slipSystemLattice(s1,i), myStructure)], [3,3])) enddo enddo @@ -1058,9 +1068,9 @@ if (.not. phase_localConstitution(phase)) then if (neighboring_ns == ns) then if (neighboring_el /= el .or. neighboring_ip /= ip) then connection_latticeConf(1:3,n) = math_mul33x3(invFe, neighboring_ipCoords - ipCoords) - forall (s = 1:ns, c = 1:2) & - neighboring_rhoExcess(c,s,n) = state(g,neighboring_ip,neighboring_el)%p((2*c-2)*ns+s) & ! positive mobiles - - state(g,neighboring_ip,neighboring_el)%p((2*c-1)*ns+s) ! negative mobiles + forall (s = 1_pInt:ns, c = 1_pInt:2_pInt) & + neighboring_rhoExcess(c,s,n) = state(g,neighboring_ip,neighboring_el)%p((2_pInt*c-2_pInt)*ns+s) & ! positive mobiles + - state(g,neighboring_ip,neighboring_el)%p((2_pInt*c-1_pInt)*ns+s) ! negative mobiles else ! thats myself! probably using periodic images connection_latticeConf(1:3,n) = areaNormal_latticeConf(1:3,n) * FVsize @@ -1068,7 +1078,7 @@ if (.not. phase_localConstitution(phase)) then endif else ! different number of active slip systems - call IO_error(-1,ext_msg='different number of active slip systems in neighboring IPs of same crystal structure') + call IO_error(-1_pInt,ext_msg='different number of active slip systems in neighboring IPs of same crystal structure') endif else ! local neighbor or different lattice structure or different constitution instance @@ -1449,15 +1459,15 @@ enddo !*** get dislocation velocity and its tangent and store the velocity in the state array if (myStructure == 1_pInt) then ! for fcc all velcities are equal - call constitutive_nonlocal_kinetics(v(1:ns,1), tau, 1, Temperature, state, g, ip, el, dv_dtau(1:ns,1)) - do t = 1,4 + call constitutive_nonlocal_kinetics(v(1:ns,1), tau, 1_pInt, Temperature, state, g, ip, el, dv_dtau(1:ns,1)) + do t = 1_pInt,4_pInt v(1:ns,t) = v(1:ns,1) dv_dtau(1:ns,t) = dv_dtau(1:ns,1) - state%p((12+t)*ns+1:(13+t)*ns) = v(1:ns,1) + state%p((12_pInt+t)*ns+1:(13_pInt+t)*ns) = v(1:ns,1) enddo else ! for all other lattice structures the velcities may vary with character and sign - do t = 1,4 - c = (t-1)/2+1 + do t = 1_pInt,4_pInt + c = (t-1_pInt)/2_pInt+1_pInt call constitutive_nonlocal_kinetics(v(1:ns,t), tau, c, Temperature, state, g, ip, el, dv_dtau(1:ns,t)) state%p((12+t)*ns+1:(13+t)*ns) = v(1:ns,t) enddo @@ -1466,8 +1476,8 @@ endif !*** Bauschinger effect -forall (s = 1:ns, t = 5:8, state%p((t-1)*ns+s) * v(s,t-4) < 0.0_pReal) & - rhoSgl(s,t-4) = rhoSgl(s,t-4) + abs(state%p((t-1)*ns+s)) +forall (s = 1_pInt:ns, t = 5_pInt:8_pInt, state%p((t-1)*ns+s) * v(s,t-4_pInt) < 0.0_pReal) & + rhoSgl(s,t-4_pInt) = rhoSgl(s,t-4_pInt) + abs(state%p((t-1_pInt)*ns+s)) !*** Calculation of gdot and its tangent @@ -1636,7 +1646,7 @@ logical considerEnteringFlux, & considerLeavingFlux #ifndef _OPENMP - if (debug_verbosity > 6 .and. ((debug_e == el .and. debug_i == ip .and. debug_g == g) .or. .not. debug_selectiveDebugger)) then + if (debug_verbosity > 6_pInt .and. ((debug_e == el .and. debug_i == ip .and. debug_g == g) .or. .not. debug_selectiveDebugger)) then write(6,*) write(6,'(a,i8,1x,i2,1x,i1)') '<< CONST >> nonlocal_dotState at el ip g ',el,ip,g write(6,*) @@ -1644,10 +1654,10 @@ logical considerEnteringFlux, & #endif select case(mesh_element(2,el)) - case (1,6,7,8,9) + case (1_pInt,6_pInt,7_pInt,8_pInt,9_pInt) ! all fine case default - call IO_error(-1,el,ip,g,'element type not supported for nonlocal constitution') + call IO_error(-1_pInt,el,ip,g,'element type not supported for nonlocal constitution') end select myInstance = phase_constitutionInstance(material_phase(g,ip,el)) @@ -1662,17 +1672,17 @@ dUpper = 0.0_pReal !*** shortcut to state variables -forall (s = 1:ns, t = 1:4) & - rhoSgl(s,t) = max(state(g,ip,el)%p((t-1)*ns+s), 0.0_pReal) -forall (s = 1:ns, t = 5:8) & - rhoSgl(s,t) = state(g,ip,el)%p((t-1)*ns+s) -forall (s = 1:ns, c = 1:2) & - rhoDip(s,c) = max(state(g,ip,el)%p((7+c)*ns+s), 0.0_pReal) -rhoForest = state(g,ip,el)%p(10*ns+1:11*ns) -tauThreshold = state(g,ip,el)%p(11*ns+1:12*ns) -tauBack = state(g,ip,el)%p(12*ns+1:13*ns) -forall (t = 1:4) & - v(1:ns,t) = state(g,ip,el)%p((12+t)*ns+1:(13+t)*ns) +forall (s = 1_pInt:ns, t = 1_pInt:4_pInt) & + rhoSgl(s,t) = max(state(g,ip,el)%p((t-1_pInt)*ns+s), 0.0_pReal) +forall (s = 1_pInt:ns, t = 5_pInt:8_pInt) & + rhoSgl(s,t) = state(g,ip,el)%p((t-1_pInt)*ns+s) +forall (s = 1_pInt:ns, c = 1_pInt:2_pInt) & + rhoDip(s,c) = max(state(g,ip,el)%p((7_pInt+c)*ns+s), 0.0_pReal) +rhoForest = state(g,ip,el)%p(10_pInt*ns+1:11_pInt*ns) +tauThreshold = state(g,ip,el)%p(11_pInt*ns+1_pInt:12_pInt*ns) +tauBack = state(g,ip,el)%p(12_pInt*ns+1:13_pInt*ns) +forall (t = 1_pInt:4_pInt) & + v(1_pInt:ns,t) = state(g,ip,el)%p((12_pInt+t)*ns+1_pInt:(13_pInt+t)*ns) !*** sanity check for timestep @@ -1687,13 +1697,13 @@ endif !**************************************************************************** !*** Calculate shear rate -forall (t = 1:4) & - gdot(1:ns,t) = rhoSgl(1:ns,t) * constitutive_nonlocal_burgers(1:ns,myInstance) * v(1:ns,t) -forall (s = 1:ns, t = 1:4, rhoSgl(s,t+4) * v(s,t) < 0.0_pReal) & ! contribution of used rho for changing sign of v +forall (t = 1_pInt:4_pInt) & + gdot(1_pInt:ns,t) = rhoSgl(1_pInt:ns,t) * constitutive_nonlocal_burgers(1:ns,myInstance) * v(1:ns,t) +forall (s = 1_pInt:ns, t = 1_pInt:4_pInt, rhoSgl(s,t+4_pInt) * v(s,t) < 0.0_pReal) & ! contribution of used rho for changing sign of v gdot(s,t) = gdot(s,t) + abs(rhoSgl(s,t+4)) * constitutive_nonlocal_burgers(s,myInstance) * v(s,t) #ifndef _OPENMP - if (debug_verbosity > 6 .and. ((debug_e == el .and. debug_i == ip .and. debug_g == g) .or. .not. debug_selectiveDebugger)) then + if (debug_verbosity > 6_pInt .and. ((debug_e == el .and. debug_i == ip .and. debug_g == g) .or. .not. debug_selectiveDebugger)) then write(6,'(a,/,10(12x,12(e12.5,1x),/))') '<< CONST >> rho / 1/m^2', rhoSgl, rhoDip write(6,'(a,/,4(12x,12(e12.5,1x),/))') '<< CONST >> gdot / 1/s',gdot endif @@ -1801,7 +1811,7 @@ if (.not. phase_localConstitution(material_phase(g,ip,el))) then enddo endif - opposite_n = n + mod(n,2) - mod(n+1,2) + opposite_n = n + mod(n,2_pInt) - mod(n+1_pInt,2_pInt) opposite_el = mesh_ipNeighborhood(1,opposite_n,ip,el) opposite_ip = mesh_ipNeighborhood(2,opposite_n,ip,el) diff --git a/code/constitutive_phenopowerlaw.f90 b/code/constitutive_phenopowerlaw.f90 index fce836cef..d24960aea 100644 --- a/code/constitutive_phenopowerlaw.f90 +++ b/code/constitutive_phenopowerlaw.f90 @@ -263,88 +263,93 @@ subroutine constitutive_phenopowerlaw_init(file) if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif - if (section > 0 .and. phase_constitution(section) == constitutive_phenopowerlaw_label) then ! one of my sections + if (section > 0_pInt .and. phase_constitution(section) == constitutive_phenopowerlaw_label) then ! one of my sections i = phase_constitutionInstance(section) ! which instance of my constitution is present phase positions = IO_stringPos(line,maxNchunks) - tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key select case(tag) case ('(output)') - output = output + 1 - constitutive_phenopowerlaw_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + constitutive_phenopowerlaw_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('lattice_structure') - constitutive_phenopowerlaw_structureName(i) = IO_lc(IO_stringValue(line,positions,2)) + constitutive_phenopowerlaw_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('covera_ratio') - constitutive_phenopowerlaw_CoverA(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_CoverA(i) = IO_floatValue(line,positions,2_pInt) case ('c11') - constitutive_phenopowerlaw_C11(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_C11(i) = IO_floatValue(line,positions,2_pInt) case ('c12') - constitutive_phenopowerlaw_C12(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_C12(i) = IO_floatValue(line,positions,2_pInt) case ('c13') - constitutive_phenopowerlaw_C13(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_C13(i) = IO_floatValue(line,positions,2_pInt) case ('c33') - constitutive_phenopowerlaw_C33(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_C33(i) = IO_floatValue(line,positions,2_pInt) case ('c44') - constitutive_phenopowerlaw_C44(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_C44(i) = IO_floatValue(line,positions,2_pInt) case ('nslip') - forall (j = 1:lattice_maxNslipFamily) constitutive_phenopowerlaw_Nslip(j,i) = IO_intValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily)& + constitutive_phenopowerlaw_Nslip(j,i) = IO_intValue(line,positions,1_pInt+j) case ('gdot0_slip') - constitutive_phenopowerlaw_gdot0_slip(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_gdot0_slip(i) = IO_floatValue(line,positions,2_pInt) case ('n_slip') - constitutive_phenopowerlaw_n_slip(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_n_slip(i) = IO_floatValue(line,positions,2_pInt) case ('tau0_slip') - forall (j = 1:lattice_maxNslipFamily) constitutive_phenopowerlaw_tau0_slip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily)& + constitutive_phenopowerlaw_tau0_slip(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('tausat_slip') - forall (j = 1:lattice_maxNslipFamily) constitutive_phenopowerlaw_tausat_slip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily)& + constitutive_phenopowerlaw_tausat_slip(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('a_slip', 'w0_slip') - constitutive_phenopowerlaw_a_slip(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_a_slip(i) = IO_floatValue(line,positions,2_pInt) case ('ntwin') - forall (j = 1:lattice_maxNtwinFamily) constitutive_phenopowerlaw_Ntwin(j,i) = IO_intValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily)& + constitutive_phenopowerlaw_Ntwin(j,i) = IO_intValue(line,positions,1_pInt+j) case ('gdot0_twin') - constitutive_phenopowerlaw_gdot0_twin(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_gdot0_twin(i) = IO_floatValue(line,positions,2_pInt) case ('n_twin') - constitutive_phenopowerlaw_n_twin(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_n_twin(i) = IO_floatValue(line,positions,2_pInt) case ('tau0_twin') - forall (j = 1:lattice_maxNtwinFamily) constitutive_phenopowerlaw_tau0_twin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily)& + constitutive_phenopowerlaw_tau0_twin(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('s_pr') - constitutive_phenopowerlaw_spr(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_spr(i) = IO_floatValue(line,positions,2_pInt) case ('twin_b') - constitutive_phenopowerlaw_twinB(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_twinB(i) = IO_floatValue(line,positions,2_pInt) case ('twin_c') - constitutive_phenopowerlaw_twinC(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_twinC(i) = IO_floatValue(line,positions,2_pInt) case ('twin_d') - constitutive_phenopowerlaw_twinD(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_twinD(i) = IO_floatValue(line,positions,2_pInt) case ('twin_e') - constitutive_phenopowerlaw_twinE(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_twinE(i) = IO_floatValue(line,positions,2_pInt) case ('h0_slipslip') - constitutive_phenopowerlaw_h0_slipslip(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_h0_slipslip(i) = IO_floatValue(line,positions,2_pInt) case ('h0_sliptwin') - constitutive_phenopowerlaw_h0_sliptwin(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_h0_sliptwin(i) = IO_floatValue(line,positions,2_pInt) case ('h0_twinslip') - constitutive_phenopowerlaw_h0_twinslip(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_h0_twinslip(i) = IO_floatValue(line,positions,2_pInt) case ('h0_twintwin') - constitutive_phenopowerlaw_h0_twintwin(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_h0_twintwin(i) = IO_floatValue(line,positions,2_pInt) case ('atol_resistance') - constitutive_phenopowerlaw_aTolResistance(i) = IO_floatValue(line,positions,2) + constitutive_phenopowerlaw_aTolResistance(i) = IO_floatValue(line,positions,2_pInt) case ('interaction_slipslip') - forall (j = 1:lattice_maxNinteraction) & - constitutive_phenopowerlaw_interaction_slipslip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_phenopowerlaw_interaction_slipslip(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('interaction_sliptwin') - forall (j = 1:lattice_maxNinteraction) & - constitutive_phenopowerlaw_interaction_sliptwin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_phenopowerlaw_interaction_sliptwin(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('interaction_twinslip') - forall (j = 1:lattice_maxNinteraction) & - constitutive_phenopowerlaw_interaction_twinslip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_phenopowerlaw_interaction_twinslip(j,i) = IO_floatValue(line,positions,1_pInt+j) case ('interaction_twintwin') - forall (j = 1:lattice_maxNinteraction) & - constitutive_phenopowerlaw_interaction_twintwin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_phenopowerlaw_interaction_twintwin(j,i) = IO_floatValue(line,positions,1_pInt+j) end select endif enddo -100 do i = 1,maxNinstance +100 do i = 1_pInt,maxNinstance constitutive_phenopowerlaw_structure(i) = lattice_initializeStructure(constitutive_phenopowerlaw_structureName(i), & ! get structure constitutive_phenopowerlaw_CoverA(i)) @@ -357,21 +362,21 @@ subroutine constitutive_phenopowerlaw_init(file) constitutive_phenopowerlaw_totalNslip(i) = sum(constitutive_phenopowerlaw_Nslip(:,i)) ! how many slip systems altogether constitutive_phenopowerlaw_totalNtwin(i) = sum(constitutive_phenopowerlaw_Ntwin(:,i)) ! how many twin systems altogether - if (constitutive_phenopowerlaw_structure(i) < 1 ) call IO_error(205,i) + if (constitutive_phenopowerlaw_structure(i) < 1 ) call IO_error(205_pInt,i) if (any(constitutive_phenopowerlaw_tau0_slip(:,i) < 0.0_pReal .and. & - constitutive_phenopowerlaw_Nslip(:,i) > 0)) call IO_error(210,i) - if (constitutive_phenopowerlaw_gdot0_slip(i) <= 0.0_pReal) call IO_error(211,i) - if (constitutive_phenopowerlaw_n_slip(i) <= 0.0_pReal) call IO_error(212,i) + constitutive_phenopowerlaw_Nslip(:,i) > 0)) call IO_error(210_pInt,i) + if (constitutive_phenopowerlaw_gdot0_slip(i) <= 0.0_pReal) call IO_error(211_pInt,i) + if (constitutive_phenopowerlaw_n_slip(i) <= 0.0_pReal) call IO_error(212_pInt,i) if (any(constitutive_phenopowerlaw_tausat_slip(:,i) <= 0.0_pReal .and. & - constitutive_phenopowerlaw_Nslip(:,i) > 0)) call IO_error(213,i) + constitutive_phenopowerlaw_Nslip(:,i) > 0)) call IO_error(213_pInt,i) if (any(constitutive_phenopowerlaw_a_slip(i) == 0.0_pReal .and. & - constitutive_phenopowerlaw_Nslip(:,i) > 0)) call IO_error(214,i) + constitutive_phenopowerlaw_Nslip(:,i) > 0)) call IO_error(214_pInt,i) if (any(constitutive_phenopowerlaw_tau0_twin(:,i) < 0.0_pReal .and. & - constitutive_phenopowerlaw_Ntwin(:,i) > 0)) call IO_error(210,i) + constitutive_phenopowerlaw_Ntwin(:,i) > 0)) call IO_error(210_pInt,i) if ( constitutive_phenopowerlaw_gdot0_twin(i) <= 0.0_pReal .and. & - any(constitutive_phenopowerlaw_Ntwin(:,i) > 0)) call IO_error(211,i) + any(constitutive_phenopowerlaw_Ntwin(:,i) > 0)) call IO_error(211_pInt,i) if ( constitutive_phenopowerlaw_n_twin(i) <= 0.0_pReal .and. & - any(constitutive_phenopowerlaw_Ntwin(:,i) > 0)) call IO_error(212,i) + any(constitutive_phenopowerlaw_Ntwin(:,i) > 0)) call IO_error(212_pInt,i) if (constitutive_phenopowerlaw_aTolResistance(i) <= 0.0_pReal) & constitutive_phenopowerlaw_aTolResistance(i) = 1.0_pReal ! default absolute tolerance 1 Pa @@ -395,8 +400,8 @@ subroutine constitutive_phenopowerlaw_init(file) constitutive_phenopowerlaw_hardeningMatrix_twintwin = 0.0_pReal - do i = 1,maxNinstance - do j = 1,maxval(phase_Noutput) + do i = 1_pInt,maxNinstance + do j = 1_pInt,maxval(phase_Noutput) select case(constitutive_phenopowerlaw_output(j,i)) case('resistance_slip', & 'shearrate_slip', & @@ -424,21 +429,21 @@ subroutine constitutive_phenopowerlaw_init(file) enddo constitutive_phenopowerlaw_sizeDotState(i) = constitutive_phenopowerlaw_totalNslip(i)+ & - constitutive_phenopowerlaw_totalNtwin(i)+ 2 ! s_slip, s_twin, sum(gamma), sum(f) + constitutive_phenopowerlaw_totalNtwin(i)+ 2_pInt ! s_slip, s_twin, sum(gamma), sum(f) constitutive_phenopowerlaw_sizeState(i) = constitutive_phenopowerlaw_totalNslip(i)+ & - constitutive_phenopowerlaw_totalNtwin(i)+ 2 ! s_slip, s_twin, sum(gamma), sum(f) + constitutive_phenopowerlaw_totalNtwin(i)+ 2_pInt ! s_slip, s_twin, sum(gamma), sum(f) myStructure = constitutive_phenopowerlaw_structure(i) select case (lattice_symmetryType(myStructure)) ! assign elasticity tensor - case(1) ! cubic(s) - forall(k=1:3) - forall(j=1:3) & + case(1_pInt) ! cubic(s) + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) & 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) + constitutive_phenopowerlaw_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_phenopowerlaw_C44(i) end forall - case(2) ! hex + case(2_pInt) ! 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) @@ -456,12 +461,12 @@ subroutine constitutive_phenopowerlaw_init(file) constitutive_phenopowerlaw_Cslip_66(:,:,i) = & math_Mandel3333to66(math_Voigt66to3333(constitutive_phenopowerlaw_Cslip_66(:,:,i))) - do f = 1,lattice_maxNslipFamily ! >>> interaction slip -- X - index_myFamily = sum(constitutive_phenopowerlaw_Nslip(1:f-1,i)) - do j = 1,constitutive_phenopowerlaw_Nslip(f,i) ! loop over (active) systems in my family (slip) - do o = 1,lattice_maxNslipFamily - index_otherFamily = sum(constitutive_phenopowerlaw_Nslip(1:o-1,i)) - do k = 1,constitutive_phenopowerlaw_Nslip(o,i) ! loop over (active) systems in other family (slip) + do f = 1_pInt,lattice_maxNslipFamily ! >>> interaction slip -- X + index_myFamily = sum(constitutive_phenopowerlaw_Nslip(1:f-1_pInt,i)) + do j = 1_pInt,constitutive_phenopowerlaw_Nslip(f,i) ! loop over (active) systems in my family (slip) + do o = 1_pInt,lattice_maxNslipFamily + index_otherFamily = sum(constitutive_phenopowerlaw_Nslip(1:o-1_pInt,i)) + do k = 1_pInt,constitutive_phenopowerlaw_Nslip(o,i) ! loop over (active) systems in other family (slip) constitutive_phenopowerlaw_hardeningMatrix_slipslip(index_otherFamily+k,index_myFamily+j,i) = & constitutive_phenopowerlaw_interaction_slipslip(lattice_interactionSlipSlip( & sum(lattice_NslipSystem(1:o-1,myStructure))+k, & @@ -469,39 +474,39 @@ subroutine constitutive_phenopowerlaw_init(file) myStructure), i ) enddo; enddo - do o = 1,lattice_maxNtwinFamily - index_otherFamily = sum(constitutive_phenopowerlaw_Ntwin(1:o-1,i)) - do k = 1,constitutive_phenopowerlaw_Ntwin(o,i) ! loop over (active) systems in other family (twin) + do o = 1_pInt,lattice_maxNtwinFamily + index_otherFamily = sum(constitutive_phenopowerlaw_Ntwin(1:o-1_pInt,i)) + do k = 1_pInt,constitutive_phenopowerlaw_Ntwin(o,i) ! loop over (active) systems in other family (twin) constitutive_phenopowerlaw_hardeningMatrix_sliptwin(index_otherFamily+k,index_myFamily+j,i) = & constitutive_phenopowerlaw_interaction_sliptwin(lattice_interactionSlipTwin( & - sum(lattice_NtwinSystem(1:o-1,myStructure))+k, & - sum(lattice_NslipSystem(1:f-1,myStructure))+j, & + sum(lattice_NtwinSystem(1:o-1_pInt,myStructure))+k, & + sum(lattice_NslipSystem(1:f-1_pInt,myStructure))+j, & myStructure), i ) enddo; enddo enddo; enddo - do f = 1,lattice_maxNtwinFamily ! >>> interaction twin -- X - index_myFamily = sum(constitutive_phenopowerlaw_Ntwin(1:f-1,i)) - do j = 1,constitutive_phenopowerlaw_Ntwin(f,i) ! loop over (active) systems in my family (twin) + do f = 1_pInt,lattice_maxNtwinFamily ! >>> interaction twin -- X + index_myFamily = sum(constitutive_phenopowerlaw_Ntwin(1:f-1_pInt,i)) + do j = 1_pInt,constitutive_phenopowerlaw_Ntwin(f,i) ! loop over (active) systems in my family (twin) - do o = 1,lattice_maxNslipFamily - index_otherFamily = sum(constitutive_phenopowerlaw_Nslip(1:o-1,i)) - do k = 1,constitutive_phenopowerlaw_Nslip(o,i) ! loop over (active) systems in other family (slip) + do o = 1_pInt,lattice_maxNslipFamily + index_otherFamily = sum(constitutive_phenopowerlaw_Nslip(1:o-1_pInt,i)) + do k = 1_pInt,constitutive_phenopowerlaw_Nslip(o,i) ! loop over (active) systems in other family (slip) constitutive_phenopowerlaw_hardeningMatrix_twinslip(index_otherFamily+k,index_myFamily+j,i) = & constitutive_phenopowerlaw_interaction_twinslip(lattice_interactionTwinSlip( & - sum(lattice_NslipSystem(1:o-1,myStructure))+k, & - sum(lattice_NtwinSystem(1:f-1,myStructure))+j, & + sum(lattice_NslipSystem(1:o-1_pInt,myStructure))+k, & + sum(lattice_NtwinSystem(1:f-1_pInt,myStructure))+j, & myStructure), i ) enddo; enddo - do o = 1,lattice_maxNtwinFamily - index_otherFamily = sum(constitutive_phenopowerlaw_Ntwin(1:o-1,i)) - do k = 1,constitutive_phenopowerlaw_Ntwin(o,i) ! loop over (active) systems in other family (twin) + do o = 1_pInt,lattice_maxNtwinFamily + index_otherFamily = sum(constitutive_phenopowerlaw_Ntwin(1:o-1_pInt,i)) + do k = 1_pInt,constitutive_phenopowerlaw_Ntwin(o,i) ! loop over (active) systems in other family (twin) constitutive_phenopowerlaw_hardeningMatrix_twintwin(index_otherFamily+k,index_myFamily+j,i) = & constitutive_phenopowerlaw_interaction_twintwin(lattice_interactionTwinTwin( & - sum(lattice_NtwinSystem(1:o-1,myStructure))+k, & - sum(lattice_NtwinSystem(1:f-1,myStructure))+j, & + sum(lattice_NtwinSystem(1:o-1_pInt,myStructure))+k, & + sum(lattice_NtwinSystem(1:f-1_pInt,myStructure))+j, & myStructure), i ) enddo; enddo diff --git a/code/constitutive_titanmod.f90 b/code/constitutive_titanmod.f90 index fff006ae3..3e2ca49e7 100644 --- a/code/constitutive_titanmod.f90 +++ b/code/constitutive_titanmod.f90 @@ -389,199 +389,199 @@ do ! read thru sections of if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif - if (section > 0 .and. phase_constitution(section) == constitutive_titanmod_label) then ! one of my sections + if (section > 0_pInt .and. phase_constitution(section) == constitutive_titanmod_label) then ! one of my sections i = phase_constitutionInstance(section) ! which instance of my constitution is present phase positions = IO_stringPos(line,maxNchunks) - tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key select case(tag) case ('(output)') - output = output + 1 - constitutive_titanmod_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + constitutive_titanmod_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) write(6,*) tag,constitutive_titanmod_output(output,i) case ('lattice_structure') - constitutive_titanmod_structureName(i) = IO_lc(IO_stringValue(line,positions,2)) + constitutive_titanmod_structureName(i) = IO_lc(IO_stringValue(line,positions,2_pInt)) write(6,*) tag,constitutive_titanmod_structureName(i) case ('covera_ratio') - constitutive_titanmod_CoverA(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_CoverA(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_CoverA(i) case ('c11') - constitutive_titanmod_C11(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_C11(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_C11(i) case ('c12') - constitutive_titanmod_C12(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_C12(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_C12(i) case ('c13') - constitutive_titanmod_C13(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_C13(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_C13(i) case ('c33') - constitutive_titanmod_C33(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_C33(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_C33(i) case ('c44') - constitutive_titanmod_C44(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_C44(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_C44(i) case ('debyefrequency') - constitutive_titanmod_debyefrequency(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_debyefrequency(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_debyefrequency(i) case ('kinkf0') - constitutive_titanmod_kinkf0(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_kinkf0(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag,constitutive_titanmod_kinkf0(i) case ('nslip') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_Nslip(j,i) = IO_intValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_Nslip(j,i) = IO_intValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_Nslip(1:4,i) case ('ntwin') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_Ntwin(j,i) = IO_intValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_Ntwin(j,i) = IO_intValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_Ntwin(1:4,i) case ('rho_edge0') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_rho_edge0(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_rho_edge0(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_rho_edge0(1:4,i) case ('rho_screw0') forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_rho_screw0(j,i) = IO_floatValue(line,positions,1+j) + constitutive_titanmod_rho_screw0(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_rho_screw0(1:4,i) case ('slipburgers') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_burgersPerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_burgersPerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_burgersPerSlipFamily(1:4,i) case ('twinburgers') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_burgersPerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_burgersPerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_burgersPerTwinFamily(1:4,i) case ('f0') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_f0_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_f0_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_f0_PerSlipFamily(1:4,i) case ('twinf0') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twinf0_PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twinf0_PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twinf0_PerTwinFamily(1:4,i) case ('tau0e') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_tau0e_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_tau0e_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_tau0e_PerSlipFamily(1:4,i) case ('twintau0') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twintau0_PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twintau0_PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twintau0_PerTwinFamily(1:4,i) case ('tau0s') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_tau0s_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_tau0s_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_tau0s_PerSlipFamily(1:4,i) case ('capre') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_capre_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_capre_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_capre_PerSlipFamily(1:4,i) case ('caprs') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_caprs_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_caprs_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_caprs_PerSlipFamily(1:4,i) case ('v0e') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_v0e_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_v0e_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_v0e_PerSlipFamily(1:4,i) case ('twingamma0') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twingamma0_PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twingamma0_PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twingamma0_PerTwinFamily(1:4,i) case ('v0s') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_v0s_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_v0s_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_v0s_PerSlipFamily(1:4,i) case ('kinkcriticallength') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_kinkcriticallength_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_kinkcriticallength_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_kinkcriticallength_PerSlipFamily(1:4,i) case ('twinsize') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twinsizePerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twinsizePerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('celambdaslip') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_CeLambdaSlipPerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_CeLambdaSlipPerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('twinlambdaslip') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twinlambdaslipPerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twinlambdaslipPerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twinlambdaslipPerTwinFamily(1:4,i) case ('cslambdaslip') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_CsLambdaSlipPerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_CsLambdaSlipPerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('grainsize') - constitutive_titanmod_GrainSize(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_GrainSize(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag case ('maxtwinfraction') - constitutive_titanmod_MaxTwinFraction(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_MaxTwinFraction(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag case ('pe') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_pe_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_pe_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_pe_PerSlipFamily(1:4,i) case ('twinp') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twinp_PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twinp_PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twinp_PerTwinFamily(1:4,i) case ('ps') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_ps_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_ps_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_ps_PerSlipFamily(1:4,i) case ('qe') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_qe_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_qe_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_qe_PerSlipFamily(1:4,i) case ('twinq') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twinq_PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twinq_PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twinq_PerTwinFamily(1:4,i) case ('qs') - forall (j = 1:lattice_maxNslipFamily) & - constitutive_titanmod_qs_PerSlipFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNslipFamily) & + constitutive_titanmod_qs_PerSlipFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_qs_PerSlipFamily(1:4,i) case ('twinshearconstant') - forall (j = 1:lattice_maxNtwinFamily) & - constitutive_titanmod_twinshearconstant_PerTwinFamily(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNtwinFamily) & + constitutive_titanmod_twinshearconstant_PerTwinFamily(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag,constitutive_titanmod_twinshearconstant_PerTwinFamily(1:4,i) case ('dc') - constitutive_titanmod_dc(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_dc(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag case ('twinhpconstant') - constitutive_titanmod_twinhpconstant(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_twinhpconstant(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag case ('atol_rho') - constitutive_titanmod_aTolRho(i) = IO_floatValue(line,positions,2) + constitutive_titanmod_aTolRho(i) = IO_floatValue(line,positions,2_pInt) write(6,*) tag case ('interactionslipslip') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interactionSlipSlip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interactionSlipSlip(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('interactionee') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interaction_ee(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interaction_ee(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('interactionss') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interaction_ss(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interaction_ss(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('interactiones') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interaction_es(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interaction_es(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('interactionsliptwin') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interactionSlipTwin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interactionSlipTwin(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('interactiontwinslip') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interactionTwinSlip(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interactionTwinSlip(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag case ('interactiontwintwin') - forall (j = 1:lattice_maxNinteraction) & - constitutive_titanmod_interactionTwinTwin(j,i) = IO_floatValue(line,positions,1+j) + forall (j = 1_pInt:lattice_maxNinteraction) & + constitutive_titanmod_interactionTwinTwin(j,i) = IO_floatValue(line,positions,1_pInt+j) write(6,*) tag end select endif @@ -589,43 +589,43 @@ enddo write(6,*) 'Material Property reading done' -100 do i = 1,maxNinstance +100 do i = 1_pInt,maxNinstance constitutive_titanmod_structure(i) = & lattice_initializeStructure(constitutive_titanmod_structureName(i),constitutive_titanmod_CoverA(i)) myStructure = constitutive_titanmod_structure(i) !* Sanity checks - if (myStructure < 1 .or. myStructure > 3) call IO_error(205) - if (sum(constitutive_titanmod_Nslip(:,i)) <= 0_pInt) call IO_error(207) - if (sum(constitutive_titanmod_Ntwin(:,i)) < 0_pInt) call IO_error(208) !*** + if (myStructure < 1_pInt .or. myStructure > 3_pInt) call IO_error(205_pInt) + if (sum(constitutive_titanmod_Nslip(:,i)) <= 0_pInt) call IO_error(207_pInt) + if (sum(constitutive_titanmod_Ntwin(:,i)) < 0_pInt) call IO_error(208_pInt) !*** do f = 1,lattice_maxNslipFamily if (constitutive_titanmod_Nslip(f,i) > 0_pInt) then - if (constitutive_titanmod_rho_edge0(f,i) < 0.0_pReal) call IO_error(209) - if (constitutive_titanmod_rho_screw0(f,i) < 0.0_pReal) call IO_error(210) - if (constitutive_titanmod_burgersPerSlipFamily(f,i) <= 0.0_pReal) call IO_error(211) - if (constitutive_titanmod_f0_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(212) - if (constitutive_titanmod_tau0e_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(229) - if (constitutive_titanmod_tau0s_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(233) - if (constitutive_titanmod_capre_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(234) - if (constitutive_titanmod_caprs_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(235) - if (constitutive_titanmod_v0e_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(226) - if (constitutive_titanmod_v0s_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(226) - if (constitutive_titanmod_kinkcriticallength_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(238) + if (constitutive_titanmod_rho_edge0(f,i) < 0.0_pReal) call IO_error(209_pInt) + if (constitutive_titanmod_rho_screw0(f,i) < 0.0_pReal) call IO_error(210_pInt) + if (constitutive_titanmod_burgersPerSlipFamily(f,i) <= 0.0_pReal) call IO_error(211_pInt) + if (constitutive_titanmod_f0_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(212_pInt) + if (constitutive_titanmod_tau0e_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(229_pInt) + if (constitutive_titanmod_tau0s_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(233_pInt) + if (constitutive_titanmod_capre_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(234_pInt) + if (constitutive_titanmod_caprs_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(235_pInt) + if (constitutive_titanmod_v0e_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(226_pInt) + if (constitutive_titanmod_v0s_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(226_pInt) + if (constitutive_titanmod_kinkcriticallength_PerSlipFamily(f,i) <= 0.0_pReal) call IO_error(238_pInt) endif enddo do f = 1,lattice_maxNtwinFamily if (constitutive_titanmod_Ntwin(f,i) > 0_pInt) then - if (constitutive_titanmod_burgersPerTwinFamily(f,i) <= 0.0_pReal) call IO_error(221) !*** - if (constitutive_titanmod_twinf0_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(228) - if (constitutive_titanmod_twinshearconstant_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(228) - if (constitutive_titanmod_twintau0_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(229) - if (constitutive_titanmod_twingamma0_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(226) + if (constitutive_titanmod_burgersPerTwinFamily(f,i) <= 0.0_pReal) call IO_error(221_pInt) !*** + if (constitutive_titanmod_twinf0_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(228_pInt) + if (constitutive_titanmod_twinshearconstant_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(228_pInt) + if (constitutive_titanmod_twintau0_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(229_pInt) + if (constitutive_titanmod_twingamma0_PerTwinFamily(f,i) <= 0.0_pReal) call IO_error(226_pInt) endif enddo ! if (any(constitutive_titanmod_interactionSlipSlip(1:maxval(lattice_interactionSlipSlip(:,:,myStructure)),i) < 1.0_pReal)) call IO_error(229) - if (constitutive_titanmod_dc(i) <= 0.0_pReal) call IO_error(231) - if (constitutive_titanmod_twinhpconstant(i) <= 0.0_pReal) call IO_error(232) - if (constitutive_titanmod_aTolRho(i) <= 0.0_pReal) call IO_error(233) + if (constitutive_titanmod_dc(i) <= 0.0_pReal) call IO_error(231_pInt) + if (constitutive_titanmod_twinhpconstant(i) <= 0.0_pReal) call IO_error(232_pInt) + if (constitutive_titanmod_aTolRho(i) <= 0.0_pReal) call IO_error(233_pInt) !* Determine total number of active slip or twin systems constitutive_titanmod_Nslip(:,i) = min(lattice_NslipSystem(:,myStructure),constitutive_titanmod_Nslip(:,i)) @@ -724,25 +724,25 @@ constitutive_titanmod_Ctwin_66 = 0.0_pReal constitutive_titanmod_Ctwin_3333 = 0.0_pReal write(6,*) 'Allocated slip system variables' -do i = 1,maxNinstance +do i = 1_pInt,maxNinstance myStructure = constitutive_titanmod_structure(i) !* Inverse lookup of slip system family l = 0_pInt - do f = 1,lattice_maxNslipFamily - do k = 1,constitutive_titanmod_Nslip(f,i) - l = l + 1 + do f = 1_pInt,lattice_maxNslipFamily + do k = 1_pInt,constitutive_titanmod_Nslip(f,i) + l = l + 1_pInt constitutive_titanmod_slipFamily(l,i) = f - constitutive_titanmod_slipSystemLattice(l,i) = sum(lattice_NslipSystem(1:f-1,myStructure)) + k + constitutive_titanmod_slipSystemLattice(l,i) = sum(lattice_NslipSystem(1:f-1_pInt,myStructure)) + k enddo; enddo !* Inverse lookup of twin system family l = 0_pInt - do f = 1,lattice_maxNtwinFamily - do k = 1,constitutive_titanmod_Ntwin(f,i) - l = l + 1 + do f = 1_pInt,lattice_maxNtwinFamily + do k = 1_pInt,constitutive_titanmod_Ntwin(f,i) + l = l + 1_pInt constitutive_titanmod_twinFamily(l,i) = f - constitutive_titanmod_twinSystemLattice(l,i) = sum(lattice_NtwinSystem(1:f-1,myStructure)) + k + constitutive_titanmod_twinSystemLattice(l,i) = sum(lattice_NtwinSystem(1:f-1_pInt,myStructure)) + k enddo; enddo !* Determine size of state array @@ -809,14 +809,14 @@ write(6,*) 'Determining elasticity matrix' !* Elasticity matrix and shear modulus according to material.config select case (myStructure) - case(1:2) ! cubic(s) - forall(k=1:3) - forall(j=1:3) & + case(1_pInt:2_pInt) ! cubic(s) + forall(k=1_pInt:3_pInt) + forall(j=1_pInt:3_pInt) & constitutive_titanmod_Cslip_66(k,j,i) = constitutive_titanmod_C12(i) constitutive_titanmod_Cslip_66(k,k,i) = constitutive_titanmod_C11(i) - constitutive_titanmod_Cslip_66(k+3,k+3,i) = constitutive_titanmod_C44(i) + constitutive_titanmod_Cslip_66(k+3_pInt,k+3_pInt,i) = constitutive_titanmod_C44(i) end forall - case(3:) ! all hex + case(3_pInt:) ! all hex constitutive_titanmod_Cslip_66(1,1,i) = constitutive_titanmod_C11(i) constitutive_titanmod_Cslip_66(2,2,i) = constitutive_titanmod_C11(i) constitutive_titanmod_Cslip_66(3,3,i) = constitutive_titanmod_C33(i) @@ -836,11 +836,12 @@ write(6,*) 'Determining elasticity matrix' 0.2_pReal*(constitutive_titanmod_C11(i)-constitutive_titanmod_C12(i))+0.3_pReal*constitutive_titanmod_C44(i) !* Construction of the twin elasticity matrices - do j=1,lattice_maxNtwinFamily - do k=1,constitutive_titanmod_Ntwin(j,i) - do l=1,3 ; do m=1,3 ; do n=1,3 ; do o=1,3 ; do p=1,3 ; do q=1,3 ; do r=1,3 ; do s=1,3 - constitutive_titanmod_Ctwin_3333(l,m,n,o,sum(constitutive_titanmod_Nslip(1:j-1,i))+k,i) = & - constitutive_titanmod_Ctwin_3333(l,m,n,o,sum(constitutive_titanmod_Nslip(1:j-1,i))+k,i) + & + do j=1_pInt,lattice_maxNtwinFamily + do k=1_pInt,constitutive_titanmod_Ntwin(j,i) + do l=1_pInt,3_pInt ; do m=1_pInt,3_pInt ; do n=1_pInt,3_pInt ; do o=1_pInt,3_pInt + do p=1_pInt,3_pInt ; do q=1_pInt,3_pInt ; do r=1_pInt,3_pInt ; do s=1_pInt,3_pInt + constitutive_titanmod_Ctwin_3333(l,m,n,o,sum(constitutive_titanmod_Nslip(1:j-1,i))+k,i) = & + constitutive_titanmod_Ctwin_3333(l,m,n,o,sum(constitutive_titanmod_Nslip(1:j-1,i))+k,i) + & constitutive_titanmod_Cslip_3333(p,q,r,s,i)*& lattice_Qtwin(l,p,sum(lattice_NslipSystem(1:j-1,myStructure))+k,myStructure)* & lattice_Qtwin(m,q,sum(lattice_NslipSystem(1:j-1,myStructure))+k,myStructure)* & diff --git a/code/crystallite.f90 b/code/crystallite.f90 index 797ad13d9..35b6ecc24 100644 --- a/code/crystallite.f90 +++ b/code/crystallite.f90 @@ -103,6 +103,7 @@ CONTAINS subroutine crystallite_init(Temperature) !*** variables and functions from other modules ***! +use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran 4.6 at the moment) use prec, only: pInt, & pReal use debug, only: debug_info, & @@ -221,8 +222,7 @@ allocate(crystallite_orientation(4,gMax,iMax,eMax)); crystalli allocate(crystallite_orientation0(4,gMax,iMax,eMax)); crystallite_orientation0 = 0.0_pReal allocate(crystallite_rotation(4,gMax,iMax,eMax)); crystallite_rotation = 0.0_pReal allocate(crystallite_disorientation(4,nMax,gMax,iMax,eMax)); crystallite_disorientation = 0.0_pReal -allocate(crystallite_symmetryID(gMax,iMax,eMax)); crystallite_symmetryID = 0_pInt -allocate(crystallite_localConstitution(gMax,iMax,eMax)); crystallite_localConstitution = .true. +allocate(crystallite_symmetryID(gMax,iMax,eMax)); crystallite_symmetryID = 0_pIntallocate(crystallite_localConstitution(gMax,iMax,eMax)); crystallite_localConstitution = .true. allocate(crystallite_requested(gMax,iMax,eMax)); crystallite_requested = .false. allocate(crystallite_todo(gMax,iMax,eMax)); crystallite_todo = .false. allocate(crystallite_converged(gMax,iMax,eMax)); crystallite_converged = .true. @@ -234,10 +234,10 @@ allocate(crystallite_sizePostResult(maxval(crystallite_Noutput), & if (.not. IO_open_jobFile(file,material_localFileExt)) then ! no local material configuration present... - if (.not. IO_open_file(file,material_configFile)) call IO_error(100) ! ...and cannot open material.config file + if (.not. IO_open_file(file,material_configFile)) call IO_error(100_pInt) ! ...and cannot open material.config file endif line = '' -section = 0 +section = 0_pInt do while (IO_lc(IO_getTag(line,'<','>')) /= material_partCrystallite) ! wind forward to read(file,'(a1024)',END=100) line @@ -248,38 +248,38 @@ do ! read thru sections of if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif - if (section > 0) then + if (section > 0_pInt) then positions = IO_stringPos(line,maxNchunks) - tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key select case(tag) case ('(output)') - output = output + 1 - crystallite_output(output,section) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + crystallite_output(output,section) = IO_lc(IO_stringValue(line,positions,2_pInt)) end select endif enddo 100 close(file) -do i = 1,material_Ncrystallite ! sanity checks +do i = 1_pInt,material_Ncrystallite ! sanity checks enddo -do i = 1,material_Ncrystallite - do j = 1,crystallite_Noutput(i) +do i = 1_pInt,material_Ncrystallite + do j = 1_pInt,crystallite_Noutput(i) select case(crystallite_output(j,i)) case('phase','texture','volume') - mySize = 1 + mySize = 1_pInt case('orientation','grainrotation') ! orientation as quaternion, or deviation from initial grain orientation in axis-angle form (angle in degrees) - mySize = 4 + mySize = 4_pInt case('eulerangles') ! Bunge (3-1-3) Euler angles - mySize = 3 + mySize = 3_pInt case('defgrad','f','fe','fp','lp','e','ee','p','firstpiola','1stpiola','s','tstar','secondpiola','2ndpiola') - mySize = 9 + mySize = 9_pInt case default - mySize = 0 + mySize = 0_pInt end select if (mySize > 0_pInt) then ! any meaningful output found @@ -290,7 +290,7 @@ do i = 1,material_Ncrystallite enddo crystallite_maxSizePostResults = 0_pInt -do j = 1,material_Nmicrostructure +do j = 1_pInt,material_Nmicrostructure if (microstructure_active(j)) & crystallite_maxSizePostResults = max(crystallite_maxSizePostResults,& crystallite_sizePostResults(microstructure_crystallite(j))) @@ -298,9 +298,9 @@ enddo ! write description file for crystallite output -if(.not. IO_write_jobFile(file,'outputCrystallite')) call IO_error (50) ! problems in writing file +if(.not. IO_write_jobFile(file,'outputCrystallite')) call IO_error (50_pInt) ! problems in writing file -do p = 1,material_Ncrystallite +do p = 1_pInt,material_Ncrystallite write(file,*) write(file,'(a)') '['//trim(crystallite_name(p))//']' write(file,*) @@ -619,10 +619,10 @@ do while (any(crystallite_subStep(:,:,FEsolving_execELem(1):FEsolving_execElem(2 crystallite_subTstar0_v(1:6,g,i,e) = crystallite_Tstar_v(1:6,g,i,e) ! ...2nd PK stress !$OMP FLUSH(crystallite_subF0) elseif (formerSubStep > subStepMinCryst) then ! this crystallite just converged - if (debug_verbosity > 4) then + if (debug_verbosity > 4_pInt) then !$OMP CRITICAL (distributionCrystallite) - debug_CrystalliteLoopDistribution(min(nCryst+1,NiterationCrystallite)) = & - debug_CrystalliteLoopDistribution(min(nCryst+1,NiterationCrystallite)) + 1 + debug_CrystalliteLoopDistribution(min(nCryst+1_pInt,NiterationCrystallite)) = & + debug_CrystalliteLoopDistribution(min(nCryst+1_pInt,NiterationCrystallite)) + 1_pInt !$OMP END CRITICAL (distributionCrystallite) endif endif @@ -640,7 +640,7 @@ do while (any(crystallite_subStep(:,:,FEsolving_execELem(1):FEsolving_execElem(2 ! cant restore dotState here, since not yet calculated in first cutback after initialization !$OMP FLUSH(crystallite_invFp) #ifndef _OPENMP - if (debug_verbosity > 4 & + if (debug_verbosity > 4_pInt & .and. ((e == debug_e .and. i == debug_i .and. g == debug_g) .or. .not. debug_selectiveDebugger)) then write(6,'(a,f10.8)') '<< CRYST >> cutback step in crystallite_stressAndItsTangent with new crystallite_subStep: ',& crystallite_subStep(g,i,e) @@ -672,20 +672,20 @@ do while (any(crystallite_subStep(:,:,FEsolving_execELem(1):FEsolving_execElem(2 if (any(crystallite_todo)) then select case(numerics_integrator(numerics_integrationMode)) - case(1) + case(1_pInt) call crystallite_integrateStateFPI() - case(2) + case(2_pInt) call crystallite_integrateStateEuler() - case(3) + case(3_pInt) call crystallite_integrateStateAdaptiveEuler() - case(4) + case(4_pInt) call crystallite_integrateStateRK4() - case(5) + case(5_pInt) call crystallite_integrateStateRKCK45() endselect endif - NiterationCrystallite = NiterationCrystallite + 1 + NiterationCrystallite = NiterationCrystallite + 1_pInt enddo ! cutback loop @@ -774,12 +774,12 @@ if(updateJaco) then ! --- INITIALIZE UNPERTURBED STATE --- select case(numerics_integrator(numerics_integrationMode)) - case(1) ! Fix-point method: restore to last converged state at end of subinc, since this is probably closest to perturbed state + case(1_pInt) ! Fix-point method: restore to last converged state at end of subinc, since this is probably closest to perturbed state !$OMP PARALLEL DO PRIVATE(myNgrains,mySizeState,mySizeDotState) do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - do g = 1,myNgrains + do g = 1_pInt,myNgrains mySizeState = constitutive_sizeState(g,i,e) mySizeDotState = constitutive_sizeDotState(g,i,e) constitutive_state(g,i,e)%p(1:mySizeState) = constitutive_state_backup(g,i,e)%p(1:mySizeState) @@ -792,13 +792,13 @@ if(updateJaco) then crystallite_Fe = Fe_backup crystallite_Lp = Lp_backup crystallite_Tstar_v = Tstar_v_backup - case(2,3) ! explicit Euler methods: nothing to restore (except for F), since we are only doing a stress integration step - case(4,5) ! explicit Runge-Kutta methods: restore to start of subinc, since we are doing a full integration of state and stress + case(2_pInt,3_pInt) ! explicit Euler methods: nothing to restore (except for F), since we are only doing a stress integration step + case(4_pInt,5_pInt) ! explicit Runge-Kutta methods: restore to start of subinc, since we are doing a full integration of state and stress !$OMP PARALLEL DO PRIVATE(myNgrains,mySizeState,mySizeDotState) do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - do g = 1,myNgrains + do g = 1_pInt,myNgrains mySizeState = constitutive_sizeState(g,i,e) mySizeDotState = constitutive_sizeDotState(g,i,e) constitutive_state(g,i,e)%p(1:mySizeState) = constitutive_subState0(g,i,e)%p(1:mySizeState) @@ -822,15 +822,15 @@ if(updateJaco) then where (crystallite_todo) crystallite_converged = .false. ! start out non-converged select case(numerics_integrator(numerics_integrationMode)) - case(1) + case(1_pInt) call crystallite_integrateStateFPI() - case(2) + case(2_pInt) call crystallite_integrateStateEuler() - case(3) + case(3_pInt) call crystallite_integrateStateAdaptiveEuler() - case(4) + case(4_pInt) call crystallite_integrateStateRK4() - case(5) + case(5_pInt) call crystallite_integrateStateRKCK45() end select @@ -838,12 +838,12 @@ if(updateJaco) then do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - do g = 1,myNgrains + do g = 1_pInt,myNgrains if (crystallite_requested(g,i,e) .and. crystallite_converged(g,i,e)) then ! converged state warrants stiffness update select case(perturbation) - case(1) + case(1_pInt) dPdF_perturbation1(1:3,1:3,k,l,g,i,e) = (crystallite_P(1:3,1:3,g,i,e) - P_backup(1:3,1:3,g,i,e)) / myPert ! tangent dP_ij/dFg_kl - case(2) + case(2_pInt) dPdF_perturbation2(1:3,1:3,k,l,g,i,e) = (crystallite_P(1:3,1:3,g,i,e) - P_backup(1:3,1:3,g,i,e)) / myPert ! tangent dP_ij/dFg_kl end select endif @@ -861,14 +861,14 @@ if(updateJaco) then do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - do g = 1,myNgrains + do g = 1_pInt,myNgrains if (crystallite_requested(g,i,e) .and. convergenceFlag_backup(g,i,e)) then ! central solution converged select case(pert_method) - case(1) + case(1_pInt) crystallite_dPdF(1:3,1:3,1:3,1:3,g,i,e) = dPdF_perturbation1(1:3,1:3,1:3,1:3,g,i,e) - case(2) + case(2_pInt) crystallite_dPdF(1:3,1:3,1:3,1:3,g,i,e) = dPdF_perturbation2(1:3,1:3,1:3,1:3,g,i,e) - case(3) + case(3_pInt) crystallite_dPdF(1:3,1:3,1:3,1:3,g,i,e) = 0.5_pReal* ( dPdF_perturbation1(1:3,1:3,1:3,1:3,g,i,e) & + dPdF_perturbation2(1:3,1:3,1:3,1:3,g,i,e)) end select @@ -887,7 +887,7 @@ if(updateJaco) then do e = FEsolving_execElem(1),FEsolving_execElem(2) myNgrains = homogenization_Ngrains(mesh_element(3,e)) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) - do g = 1,myNgrains + do g = 1_pInt,myNgrains mySizeState = constitutive_sizeState(g,i,e) mySizeDotState = constitutive_sizeDotState(g,i,e) constitutive_state(g,i,e)%p(1:mySizeState) = constitutive_state_backup(g,i,e)%p(1:mySizeState) @@ -977,7 +977,7 @@ else eIter = FEsolving_execElem(1:2) do e = eIter(1),eIter(2) iIter(1:2,e) = FEsolving_execIP(1:2,e) - gIter(1:2,e) = (/1,homogenization_Ngrains(mesh_element(3,e))/) + gIter(1:2,e) = [ 1_pInt,homogenization_Ngrains(mesh_element(3,e))] enddo singleRun = .false. endif @@ -1724,7 +1724,7 @@ else eIter = FEsolving_execElem(1:2) do e = eIter(1),eIter(2) iIter(1:2,e) = FEsolving_execIP(1:2,e) - gIter(1:2,e) = (/1,homogenization_Ngrains(mesh_element(3,e))/) + gIter(1:2,e) = [1_pInt,homogenization_Ngrains(mesh_element(3,e))] enddo singleRun = .false. endif @@ -1996,7 +1996,7 @@ else eIter = FEsolving_execElem(1:2) do e = eIter(1),eIter(2) iIter(1:2,e) = FEsolving_execIP(1:2,e) - gIter(1:2,e) = (/1,homogenization_Ngrains(mesh_element(3,e))/) + gIter(1:2,e) = [1_pInt,homogenization_Ngrains(mesh_element(3,e))] enddo singleRun = .false. endif @@ -2890,14 +2890,14 @@ LpLoop: do dT_dLp(3*(h-1)+j,3*(k-1)+l) = dT_dLp(3*(h-1)+j,3*(k-1)+l) + C(h,j,l,m) * AB(k,m) + C(h,j,m,l) * BTA(m,k) enddo; enddo; enddo; enddo; enddo dT_dLp = -0.5_pReal * dt * dT_dLp - dR_dLp = math_identity2nd(9) - math_mul99x99(dLp_dT_constitutive, dT_dLp) + dR_dLp = math_identity2nd(9_pInt) - math_mul99x99(dLp_dT_constitutive, dT_dLp) inv_dR_dLp = 0.0_pReal - call math_invert(9,dR_dLp,inv_dR_dLp,dummy,error) ! invert dR/dLp --> dLp/dR + call math_invert(9_pInt,dR_dLp,inv_dR_dLp,dummy,error) ! invert dR/dLp --> dLp/dR if (error) then #ifndef _OPENMP - if (debug_verbosity > 4) then + if (debug_verbosity > 4_pInt) then write(6,'(a,i8,1x,i2,1x,i3,a,i3)') '<< CRYST >> integrateStress failed on dR/dLp inversion at el ip g ',e,i,g - if (debug_verbosity > 5 & + if (debug_verbosity > 5_pInt & .and. ((e == debug_e .and. i == debug_i .and. g == debug_g) .or. .not. debug_selectiveDebugger)) then write(6,*) write(6,'(a,/,9(12x,9(e15.3,1x)/))') '<< CRYST >> dR_dLp',transpose(dR_dLp) @@ -3057,8 +3057,8 @@ logical error call math_pDecomposition(crystallite_Fe(1:3,1:3,g,i,e), U, R, error) ! polar decomposition of Fe if (error) then - call IO_warning(650, e, i, g) - orientation = (/1.0_pReal, 0.0_pReal, 0.0_pReal, 0.0_pReal/) ! fake orientation + call IO_warning(650_pInt, e, i, g) + orientation = [1.0_pReal, 0.0_pReal, 0.0_pReal, 0.0_pReal] ! fake orientation else orientation = math_RtoQuaternion(transpose(R)) endif diff --git a/code/homogenization_RGC.f90 b/code/homogenization_RGC.f90 index aa3957d7f..e54fbdf05 100644 --- a/code/homogenization_RGC.f90 +++ b/code/homogenization_RGC.f90 @@ -117,44 +117,44 @@ subroutine homogenization_RGC_init(& if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif - if (section > 0 .and. homogenization_type(section) == homogenization_RGC_label) then ! one of my sections + if (section > 0_pInt .and. homogenization_type(section) == homogenization_RGC_label) then ! one of my sections i = homogenization_typeInstance(section) ! which instance of my type is present homogenization positions = IO_stringPos(line,maxNchunks) - tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key select case(tag) case ('(output)') - output = output + 1 - homogenization_RGC_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + homogenization_RGC_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('clustersize') - homogenization_RGC_Ngrains(1,i) = IO_intValue(line,positions,2) - homogenization_RGC_Ngrains(2,i) = IO_intValue(line,positions,3) - homogenization_RGC_Ngrains(3,i) = IO_intValue(line,positions,4) + homogenization_RGC_Ngrains(1,i) = IO_intValue(line,positions,2_pInt) + homogenization_RGC_Ngrains(2,i) = IO_intValue(line,positions,3_pInt) + homogenization_RGC_Ngrains(3,i) = IO_intValue(line,positions,4_pInt) case ('scalingparameter') - homogenization_RGC_xiAlpha(i) = IO_floatValue(line,positions,2) + homogenization_RGC_xiAlpha(i) = IO_floatValue(line,positions,2_pInt) case ('overproportionality') - homogenization_RGC_ciAlpha(i) = IO_floatValue(line,positions,2) + homogenization_RGC_ciAlpha(i) = IO_floatValue(line,positions,2_pInt) case ('grainsize') - homogenization_RGC_dAlpha(1,i) = IO_floatValue(line,positions,2) - homogenization_RGC_dAlpha(2,i) = IO_floatValue(line,positions,3) - homogenization_RGC_dAlpha(3,i) = IO_floatValue(line,positions,4) + homogenization_RGC_dAlpha(1,i) = IO_floatValue(line,positions,2_pInt) + homogenization_RGC_dAlpha(2,i) = IO_floatValue(line,positions,3_pInt) + homogenization_RGC_dAlpha(3,i) = IO_floatValue(line,positions,4_pInt) case ('clusterorientation') - homogenization_RGC_angles(1,i) = IO_floatValue(line,positions,2) - homogenization_RGC_angles(2,i) = IO_floatValue(line,positions,3) - homogenization_RGC_angles(3,i) = IO_floatValue(line,positions,4) + homogenization_RGC_angles(1,i) = IO_floatValue(line,positions,2_pInt) + homogenization_RGC_angles(2,i) = IO_floatValue(line,positions,3_pInt) + homogenization_RGC_angles(3,i) = IO_floatValue(line,positions,4_pInt) end select endif enddo !*** assigning cluster orientations - do e = 1,mesh_NcpElems + do e = 1_pInt,mesh_NcpElems if (homogenization_type(mesh_element(3,e)) == homogenization_RGC_label) then myInstance = homogenization_typeInstance(mesh_element(3,e)) if (all (homogenization_RGC_angles(:,myInstance) >= 399.9_pReal)) then homogenization_RGC_orientation(:,:,1,e) = math_EulerToR(math_sampleRandomOri()) - do i = 1,FE_Nips(mesh_element(2,e)) + do i = 1_pInt,FE_Nips(mesh_element(2,e)) if (microstructure_elemhomo(mesh_element(4,e))) then homogenization_RGC_orientation(:,:,i,e) = homogenization_RGC_orientation(:,:,1,e) else @@ -162,44 +162,44 @@ subroutine homogenization_RGC_init(& endif enddo else - do i = 1,FE_Nips(mesh_element(2,e)) + do i = 1_pInt,FE_Nips(mesh_element(2,e)) homogenization_RGC_orientation(:,:,i,e) = math_EulerToR(homogenization_RGC_angles(:,myInstance)*inRad) enddo endif endif enddo -100 if (debug_verbosity == 4) then +100 if (debug_verbosity == 4_pInt) then !$OMP CRITICAL (write2out) - do i = 1,maxNinstance + do i = 1_pInt,maxNinstance write(6,'(a15,1x,i4)') 'instance: ', i write(6,*) - write(6,'(a25,3(1x,i8))') 'cluster size: ',(homogenization_RGC_Ngrains(j,i),j=1,3) + write(6,'(a25,3(1x,i8))') 'cluster size: ',(homogenization_RGC_Ngrains(j,i),j=1_pInt,3_pInt) write(6,'(a25,1x,e10.3)') 'scaling parameter: ', homogenization_RGC_xiAlpha(i) write(6,'(a25,1x,e10.3)') 'over-proportionality: ', homogenization_RGC_ciAlpha(i) - write(6,'(a25,3(1x,e10.3))') 'grain size: ',(homogenization_RGC_dAlpha(j,i),j=1,3) - write(6,'(a25,3(1x,e10.3))') 'cluster orientation: ',(homogenization_RGC_angles(j,i),j=1,3) + write(6,'(a25,3(1x,e10.3))') 'grain size: ',(homogenization_RGC_dAlpha(j,i),j=1_pInt,3_pInt) + write(6,'(a25,3(1x,e10.3))') 'cluster orientation: ',(homogenization_RGC_angles(j,i),j=1_pInt,3_pInt) enddo !$OMP END CRITICAL (write2out) endif - do i = 1,maxNinstance - do j = 1,maxval(homogenization_Noutput) + do i = 1_pInt,maxNinstance + do j = 1_pInt,maxval(homogenization_Noutput) select case(homogenization_RGC_output(j,i)) case('constitutivework') - mySize = 1 + mySize = 1_pInt case('magnitudemismatch') - mySize = 3 + mySize = 3_pInt case('penaltyenergy') - mySize = 1 + mySize = 1_pInt case('volumediscrepancy') - mySize = 1 + mySize = 1_pInt case('averagerelaxrate') - mySize = 1 + mySize = 1_pInt case('maximumrelaxrate') - mySize = 1 + mySize = 1_pInt case default - mySize = 0 + mySize = 0_pInt end select if (mySize > 0_pInt) then ! any meaningful output found @@ -212,9 +212,9 @@ subroutine homogenization_RGC_init(& homogenization_RGC_sizeState(i) & - = 3*(homogenization_RGC_Ngrains(1,i)-1)*homogenization_RGC_Ngrains(2,i)*homogenization_RGC_Ngrains(3,i) & - + 3*homogenization_RGC_Ngrains(1,i)*(homogenization_RGC_Ngrains(2,i)-1)*homogenization_RGC_Ngrains(3,i) & - + 3*homogenization_RGC_Ngrains(1,i)*homogenization_RGC_Ngrains(2,i)*(homogenization_RGC_Ngrains(3,i)-1) & + = 3*(homogenization_RGC_Ngrains(1,i)-1_pInt)*homogenization_RGC_Ngrains(2,i)*homogenization_RGC_Ngrains(3,i) & + + 3*homogenization_RGC_Ngrains(1,i)*(homogenization_RGC_Ngrains(2,i)-1_pInt)*homogenization_RGC_Ngrains(3,i) & + + 3*homogenization_RGC_Ngrains(1,i)*homogenization_RGC_Ngrains(2,i)*(homogenization_RGC_Ngrains(3,i)-1_pInt) & + 8_pInt ! (1) Average constitutive work, (2-4) Overall mismatch, (5) Average penalty energy, ! (6) Volume discrepancy, (7) Avg relaxation rate component, (8) Max relaxation rate component enddo @@ -272,16 +272,16 @@ subroutine homogenization_RGC_partitionDeformation(& integer(pInt), dimension (3) :: iGrain3 integer(pInt) homID, iGrain,iFace,i,j ! - integer(pInt), parameter :: nFace = 6 + integer(pInt), parameter :: nFace = 6_pInt !* Debugging the overall deformation gradient - if (debug_verbosity == 4) then + if (debug_verbosity == 4_pInt) then !$OMP CRITICAL (write2out) write(6,'(1x,a,i3,a,i3,a)')'========== Increment: ',theInc,' Cycle: ',cycleCounter,' ==========' write(6,'(1x,a32)')'Overall deformation gradient: ' - do i = 1,3 - write(6,'(1x,3(e14.8,1x))')(avgF(i,j), j = 1,3) + do i = 1_pInt,3_pInt + write(6,'(1x,3(e14.8,1x))')(avgF(i,j), j = 1_pInt,3_pInt) enddo write(6,*)' ' call flush(6) @@ -293,21 +293,21 @@ subroutine homogenization_RGC_partitionDeformation(& F = 0.0_pReal do iGrain = 1,homogenization_Ngrains(mesh_element(3,el)) iGrain3 = homogenization_RGC_grain1to3(iGrain,homID) - do iFace = 1,nFace + do iFace = 1_pInt,nFace intFace = homogenization_RGC_getInterface(iFace,iGrain3) ! identifying 6 interfaces of each grain aVect = homogenization_RGC_relaxationVector(intFace,state,homID)! get the relaxation vectors for each interface from global relaxation vector array nVect = homogenization_RGC_interfaceNormal(intFace,ip,el) ! get the normal of each interface - forall (i=1:3,j=1:3) & + forall (i=1_pInt:3_pInt,j=1_pInt:3_pInt) & F(i,j,iGrain) = F(i,j,iGrain) + aVect(i)*nVect(j) ! calculating deformation relaxations due to interface relaxation enddo F(:,:,iGrain) = F(:,:,iGrain) + avgF(:,:) ! resulting relaxed deformation gradient !* Debugging the grain deformation gradients - if (debug_verbosity == 4) then + if (debug_verbosity == 4_pInt) then !$OMP CRITICAL (write2out) write(6,'(1x,a32,1x,i3)')'Deformation gradient of grain: ',iGrain - do i = 1,3 - write(6,'(1x,3(e14.8,1x))')(F(i,j,iGrain), j = 1,3) + do i = 1_pInt,3_pInt + write(6,'(1x,3(e14.8,1x))')(F(i,j,iGrain), j = 1_pInt,3_pInt) enddo write(6,*)' ' call flush(6) @@ -368,7 +368,7 @@ function homogenization_RGC_updateState(& real(pReal) residMax,stresMax,constitutiveWork,penaltyEnergy,volDiscrep logical error ! - integer(pInt), parameter :: nFace = 6 + integer(pInt), parameter :: nFace = 6_pInt ! real(pReal), dimension(:,:), allocatable :: tract,jmatrix,jnverse,smatrix,pmatrix,rmatrix real(pReal), dimension(:), allocatable :: resid,relax,p_relax,p_resid,drelax @@ -383,11 +383,11 @@ function homogenization_RGC_updateState(& + nGDim(1)*nGDim(2)*(nGDim(3)-1) !* Allocate the size of the global relaxation arrays/jacobian matrices depending on the size of the cluster - allocate(resid(3*nIntFaceTot)); resid = 0.0_pReal + allocate(resid(3_pInt*nIntFaceTot)); resid = 0.0_pReal allocate(tract(nIntFaceTot,3)); tract = 0.0_pReal - allocate(relax(3*nIntFaceTot)); relax = state%p(1:3*nIntFaceTot) - allocate(drelax(3*nIntFaceTot)) - drelax = state%p(1:3*nIntFaceTot) - state0%p(1:3*nIntFaceTot) + allocate(relax(3_pInt*nIntFaceTot)); relax = state%p(1:3_pInt*nIntFaceTot) + allocate(drelax(3_pInt*nIntFaceTot)) + drelax = state%p(1:3_pInt*nIntFaceTot) - state0%p(1:3_pInt*nIntFaceTot) !* Debugging the obtained state if (debug_verbosity == 4) then @@ -432,19 +432,19 @@ function homogenization_RGC_updateState(& !* Identify the left/bottom/back grain (-|N) iGr3N = faceID(2:4) ! identifying the grain ID in local coordinate system (3-dimensional index) iGrN = homogenization_RGC_grain3to1(iGr3N,homID) ! translate the local grain ID into global coordinate system (1-dimensional index) - intFaceN = homogenization_RGC_getInterface(2*faceID(1),iGr3N) + intFaceN = homogenization_RGC_getInterface(2_pInt*faceID(1),iGr3N) normN = homogenization_RGC_interfaceNormal(intFaceN,ip,el) ! get the interface normal !* Identify the right/up/front grain (+|P) iGr3P = iGr3N - iGr3P(faceID(1)) = iGr3N(faceID(1))+1 ! identifying the grain ID in local coordinate system (3-dimensional index) + iGr3P(faceID(1)) = iGr3N(faceID(1))+1_pInt ! identifying the grain ID in local coordinate system (3-dimensional index) iGrP = homogenization_RGC_grain3to1(iGr3P,homID) ! translate the local grain ID into global coordinate system (1-dimensional index) - intFaceP = homogenization_RGC_getInterface(2*faceID(1)-1,iGr3P) + intFaceP = homogenization_RGC_getInterface(2_pInt*faceID(1)-1_pInt,iGr3P) normP = homogenization_RGC_interfaceNormal(intFaceP,ip,el) ! get the interface normal !* Compute the residual of traction at the interface (in local system, 4-dimensional index) do i = 1,3 - tract(iNum,i) = sign(viscModus_RGC*(abs(drelax(i+3*(iNum-1)))/(refRelaxRate_RGC*dt))**viscPower_RGC, & + tract(iNum,i) = sign(viscModus_RGC*(abs(drelax(i+3*(iNum-1_pInt)))/(refRelaxRate_RGC*dt))**viscPower_RGC, & drelax(i+3*(iNum-1))) ! contribution from the relaxation viscosity do j = 1,3 tract(iNum,i) = tract(iNum,i) + (P(i,j,iGrP) + R(i,j,iGrP) + D(i,j,iGrP))*normP(j) & @@ -720,15 +720,15 @@ function homogenization_RGC_updateState(& !* ------------------------------------------------------------------------------------------------------------- !*** Computing the update of the state variable (relaxation vectors) using the Jacobian matrix - allocate(jnverse(3*nIntFaceTot,3*nIntFaceTot)); jnverse = 0.0_pReal - call math_invert(3*nIntFaceTot,jmatrix,jnverse,ival,error) ! Compute the inverse of the overall Jacobian matrix + allocate(jnverse(3_pInt*nIntFaceTot,3_pInt*nIntFaceTot)); jnverse = 0.0_pReal + call math_invert(3_pInt*nIntFaceTot,jmatrix,jnverse,ival,error) ! Compute the inverse of the overall Jacobian matrix !* Debugging the inverse Jacobian matrix - if (debug_verbosity == 4) then + if (debug_verbosity == 4_pInt) then !$OMP CRITICAL (write2out) write(6,'(1x,a30)')'Jacobian inverse' do i = 1,3*nIntFaceTot - write(6,'(1x,100(e10.4,1x))')(jnverse(i,j), j = 1,3*nIntFaceTot) + write(6,'(1x,100(e10.4,1x))')(jnverse(i,j), j = 1_pInt,3_pInt*nIntFaceTot) enddo write(6,*)' ' call flush(6) @@ -1292,9 +1292,9 @@ function homogenization_RGC_grain1to3(& !* Get the grain position nGDim = homogenization_RGC_Ngrains(:,homID) - homogenization_RGC_grain1to3(3) = 1+(grain1-1)/(nGDim(1)*nGDim(2)) - homogenization_RGC_grain1to3(2) = 1+mod((grain1-1)/nGDim(1),nGDim(2)) - homogenization_RGC_grain1to3(1) = 1+mod((grain1-1),nGDim(1)) + homogenization_RGC_grain1to3(3) = 1_pInt+(grain1-1_pInt)/(nGDim(1)*nGDim(2)) + homogenization_RGC_grain1to3(2) = 1_pInt+mod((grain1-1_pInt)/nGDim(1),nGDim(2)) + homogenization_RGC_grain1to3(1) = 1_pInt+mod((grain1-1_pInt),nGDim(1)) endfunction @@ -1390,21 +1390,23 @@ function homogenization_RGC_interface1to4(& !* Get the corresponding interface ID in 4D (normal and local position) if (iFace1D > 0 .and. iFace1D <= nIntFace(1)) then ! ... interface with normal //e1 homogenization_RGC_interface1to4(1) = 1 - homogenization_RGC_interface1to4(3) = mod((iFace1D-1),nGDim(2))+1 - homogenization_RGC_interface1to4(4) = mod(int(dble(iFace1D-1)/dble(nGDim(2))),nGDim(3))+1 - homogenization_RGC_interface1to4(2) = int(dble(iFace1D-1)/dble(nGDim(2))/dble(nGDim(3)))+1 + homogenization_RGC_interface1to4(3) = mod((iFace1D-1_pInt),nGDim(2))+1_pInt + homogenization_RGC_interface1to4(4) = mod(int(real(iFace1D-1_pInt,pReal)/real(nGDim(2),pReal),pInt),nGDim(3))+1 + homogenization_RGC_interface1to4(2) = int(real(iFace1D-1,pReal)/real(nGDim(2),pReal)/real(nGDim(3),pReal),pInt)+1 elseif (iFace1D > nIntFace(1) .and. iFace1D <= (nIntFace(2) + nIntFace(1))) then ! ... interface with normal //e2 homogenization_RGC_interface1to4(1) = 2 homogenization_RGC_interface1to4(4) = mod((iFace1D-nIntFace(1)-1),nGDim(3))+1 - homogenization_RGC_interface1to4(2) = mod(int(dble(iFace1D-nIntFace(1)-1)/dble(nGDim(3))),nGDim(1))+1 - homogenization_RGC_interface1to4(3) = int(dble(iFace1D-nIntFace(1)-1)/dble(nGDim(3))/dble(nGDim(1)))+1 + homogenization_RGC_interface1to4(2) = mod(int(real(iFace1D-nIntFace(1)-1,pReal)/real(nGDim(3),pReal)),nGDim(1))+1 + homogenization_RGC_interface1to4(3) = int(real(iFace1D-nIntFace(1)-1,pReal)/real(nGDim(3),pReal)/real(nGDim(1),pReal))+1 elseif (iFace1D > nIntFace(2) + nIntFace(1) .and. iFace1D <= (nIntFace(3) + nIntFace(2) + nIntFace(1))) then ! ... interface with normal //e3 homogenization_RGC_interface1to4(1) = 3 homogenization_RGC_interface1to4(2) = mod((iFace1D-nIntFace(2)-nIntFace(1)-1),nGDim(1))+1 - homogenization_RGC_interface1to4(3) = mod(int(dble(iFace1D-nIntFace(2)-nIntFace(1)-1)/dble(nGDim(1))),nGDim(2))+1 - homogenization_RGC_interface1to4(4) = int(dble(iFace1D-nIntFace(2)-nIntFace(1)-1)/dble(nGDim(1))/dble(nGDim(2)))+1 + homogenization_RGC_interface1to4(3) = mod(int(real(iFace1D-nIntFace(2)-nIntFace(1)-1,pReal)/& + real(nGDim(1),pReal),pInt),nGDim(2))+1 + homogenization_RGC_interface1to4(4) = int(real(iFace1D-nIntFace(2)-nIntFace(1)-1,pReal)/& + real(nGDim(1),pReal)/real(nGDim(2),pReal),pInt)+1 endif endfunction diff --git a/code/homogenization_isostrain.f90 b/code/homogenization_isostrain.f90 index 5512fce13..e3c82e389 100644 --- a/code/homogenization_isostrain.f90 +++ b/code/homogenization_isostrain.f90 @@ -81,7 +81,7 @@ subroutine homogenization_isostrain_init(& !$OMP END CRITICAL (write2out) maxNinstance = count(homogenization_type == homogenization_isostrain_label) - if (maxNinstance == 0) return + if (maxNinstance == 0_pInt) return allocate(homogenization_isostrain_sizeState(maxNinstance)) ; homogenization_isostrain_sizeState = 0_pInt allocate(homogenization_isostrain_sizePostResults(maxNinstance)); homogenization_isostrain_sizePostResults = 0_pInt @@ -93,7 +93,7 @@ subroutine homogenization_isostrain_init(& rewind(file) line = '' - section = 0 + section = 0_pInt do while (IO_lc(IO_getTag(line,'<','>')) /= material_partHomogenization) ! wind forward to read(file,'(a1024)',END=100) line @@ -104,24 +104,24 @@ subroutine homogenization_isostrain_init(& if (IO_isBlank(line)) cycle ! skip empty lines if (IO_getTag(line,'<','>') /= '') exit ! stop at next part if (IO_getTag(line,'[',']') /= '') then ! next section - section = section + 1 - output = 0 ! reset output counter + section = section + 1_pInt + output = 0_pInt ! reset output counter endif if (section > 0 .and. homogenization_type(section) == homogenization_isostrain_label) then ! one of my sections i = homogenization_typeInstance(section) ! which instance of my type is present homogenization positions = IO_stringPos(line,maxNchunks) - tag = IO_lc(IO_stringValue(line,positions,1)) ! extract key + tag = IO_lc(IO_stringValue(line,positions,1_pInt)) ! extract key select case(tag) case ('(output)') - output = output + 1 - homogenization_isostrain_output(output,i) = IO_lc(IO_stringValue(line,positions,2)) + output = output + 1_pInt + homogenization_isostrain_output(output,i) = IO_lc(IO_stringValue(line,positions,2_pInt)) case ('ngrains') - homogenization_isostrain_Ngrains(i) = IO_intValue(line,positions,2) + homogenization_isostrain_Ngrains(i) = IO_intValue(line,positions,2_pInt) end select endif enddo -100 do i = 1,maxNinstance ! sanity checks +100 do i = 1_pInt,maxNinstance ! sanity checks enddo do i = 1,maxNinstance diff --git a/processing/post/postResults.py b/processing/post/postResults.py index d0fba76ab..0a43460d6 100755 --- a/processing/post/postResults.py +++ b/processing/post/postResults.py @@ -93,7 +93,7 @@ class MPIEspectral_result: # mimic py_post result object self.N_loadcases = self._keyedPackedArray('loadcases',count=1,type='i',default=1)[0] self._frequencies = self._keyedPackedArray('frequencies',count=self.N_loadcases,type='i',default=1) self._increments = self._keyedPackedArray('increments',count=self.N_loadcases,type='i') - self._increments[0] -= 1 # delete zero'th entry + #self._increments[0] -= 1 # delete zero'th entry (might be needed for older spectralOuts self._times = self._keyedPackedArray('times',count=self.N_loadcases,type='d',default=0.0) self._logscales = self._keyedPackedArray('logscales',count=self.N_loadcases,type='i',default=0) self.dimension = self._keyedPackedArray('dimension',count=3,type='d')