From 194824fd0fa25a581586428c9630797c36429c85 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Fri, 22 Feb 2019 20:37:41 +0100 Subject: [PATCH] WIP: cleaned no file reading getting rid of a number of obsolete dependencies --- src/constitutive.f90 | 2 +- src/crystallite.f90 | 2 +- src/kinematics_thermal_expansion.f90 | 229 ++++++++++++----------- src/material.f90 | 8 +- src/numerics.f90 | 10 +- src/source_thermal_dissipation.f90 | 1 + src/source_thermal_externalheat.f90 | 270 ++++++++++----------------- src/thermal_adiabatic.f90 | 2 - src/thermal_conduction.f90 | 2 - 9 files changed, 225 insertions(+), 301 deletions(-) diff --git a/src/constitutive.f90 b/src/constitutive.f90 index 48d4a3f8f..a340d0adf 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -153,7 +153,7 @@ subroutine constitutive_init() ! parse source mechanisms from config file call IO_checkAndRewind(FILEUNIT) if (any(phase_source == SOURCE_thermal_dissipation_ID)) call source_thermal_dissipation_init - if (any(phase_source == SOURCE_thermal_externalheat_ID)) call source_thermal_externalheat_init(FILEUNIT) + if (any(phase_source == SOURCE_thermal_externalheat_ID)) call source_thermal_externalheat_init if (any(phase_source == SOURCE_damage_isoBrittle_ID)) call source_damage_isoBrittle_init if (any(phase_source == SOURCE_damage_isoDuctile_ID)) call source_damage_isoDuctile_init if (any(phase_source == SOURCE_damage_anisoBrittle_ID)) call source_damage_anisoBrittle_init diff --git a/src/crystallite.f90 b/src/crystallite.f90 index c272abd07..ec9c782ea 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -248,7 +248,7 @@ subroutine crystallite_init allocate(crystallite_sizePostResult(maxval(crystallite_Noutput), & size(config_crystallite)), source=0_pInt) - select case(numerics_integrator(1)) + select case(numerics_integrator) case(1_pInt) integrateState => integrateStateFPI case(2_pInt) diff --git a/src/kinematics_thermal_expansion.f90 b/src/kinematics_thermal_expansion.f90 index a44bc6902..56caa6e4b 100644 --- a/src/kinematics_thermal_expansion.f90 +++ b/src/kinematics_thermal_expansion.f90 @@ -4,21 +4,24 @@ !> @details to be done !-------------------------------------------------------------------------------------------------- module kinematics_thermal_expansion - use prec, only: & - pReal, & - pInt + use prec, only: & + pReal, & + pInt - implicit none - private + implicit none + private - !type, private :: tParameters - ! real(pReal), allocatable, dimension(:) :: & - !end type tParameters + type, private :: tParameters + real(pReal), allocatable, dimension(:,:,:) :: & + expansion + end type tParameters - public :: & - kinematics_thermal_expansion_init, & - kinematics_thermal_expansion_initialStrain, & - kinematics_thermal_expansion_LiAndItsTangent + type(tParameters), dimension(:), allocatable :: param + + public :: & + kinematics_thermal_expansion_init, & + kinematics_thermal_expansion_initialStrain, & + kinematics_thermal_expansion_LiAndItsTangent contains @@ -28,120 +31,128 @@ contains !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- subroutine kinematics_thermal_expansion_init() -#if defined(__GFORTRAN__) || __INTEL_COMPILER >= 1800 - use, intrinsic :: iso_fortran_env, only: & - compiler_version, & - compiler_options -#endif - use debug, only: & - debug_level,& - debug_constitutive,& - debug_levelBasic - use IO, only: & - IO_timeStamp - use material, only: & - phase_kinematics, & - KINEMATICS_thermal_expansion_label, & - KINEMATICS_thermal_expansion_ID - use config, only: & - config_phase - - implicit none - integer(pInt) :: & - Ninstance, & - p + use debug, only: & + debug_level,& + debug_constitutive,& + debug_levelBasic + use material, only: & + phase_kinematics, & + KINEMATICS_thermal_expansion_label, & + KINEMATICS_thermal_expansion_ID + use config, only: & + config_phase - write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_thermal_expansion_LABEL//' init -+>>>' - write(6,'(a15,a)') ' Current time: ',IO_timeStamp() -#include "compilation_info.f90" - - Ninstance = int(count(phase_kinematics == KINEMATICS_thermal_expansion_ID),pInt) - if (Ninstance == 0_pInt) return + implicit none + integer(pInt) :: & + Ninstance, & + p, i + real(pReal), dimension(:), allocatable :: & + temp + + write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_thermal_expansion_LABEL//' init -+>>>' - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & - write(6,'(a16,1x,i5,/)') '# instances:',Ninstance - - do p = 1_pInt, size(phase_kinematics) - if (all(phase_kinematics(:,p) /= KINEMATICS_thermal_expansion_ID)) cycle - enddo + Ninstance = int(count(phase_kinematics == KINEMATICS_thermal_expansion_ID),pInt) + + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & + write(6,'(a16,1x,i5,/)') '# instances:',Ninstance + + allocate(param(Ninstance)) + + do p = 1_pInt, size(phase_kinematics) + if (all(phase_kinematics(:,p) /= KINEMATICS_thermal_expansion_ID)) cycle + + ! ToDo: Here we need to decide how to extend the concept of instances to + ! kinetics and sources. I would suggest that the same mechanism exists at maximum once per phase + + ! read up to three parameters (constant, linear, quadratic with T) + temp = config_phase(p)%getFloats('thermal_expansion11') + !lattice_thermalExpansion33(1,1,1:size(temp),p) = temp + temp = config_phase(p)%getFloats('thermal_expansion22', & + defaultVal=[(0.0_pReal, i=1,size(temp))],requiredSize=size(temp)) + !lattice_thermalExpansion33(2,2,1:size(temp),p) = temp + temp = config_phase(p)%getFloats('thermal_expansion33', & + defaultVal=[(0.0_pReal, i=1,size(temp))],requiredSize=size(temp)) + enddo end subroutine kinematics_thermal_expansion_init + !-------------------------------------------------------------------------------------------------- !> @brief report initial thermal strain based on current temperature deviation from reference !-------------------------------------------------------------------------------------------------- pure function kinematics_thermal_expansion_initialStrain(homog,phase,offset) - use material, only: & - temperature - use lattice, only: & - lattice_thermalExpansion33, & - lattice_referenceTemperature + use material, only: & + temperature + use lattice, only: & + lattice_thermalExpansion33, & + lattice_referenceTemperature + + implicit none + integer(pInt), intent(in) :: & + phase, & + homog, offset + real(pReal), dimension(3,3) :: & + kinematics_thermal_expansion_initialStrain !< initial thermal strain (should be small strain, though) - implicit none - integer(pInt), intent(in) :: & - phase, & - homog, offset - real(pReal), dimension(3,3) :: & - kinematics_thermal_expansion_initialStrain !< initial thermal strain (should be small strain, though) - - - kinematics_thermal_expansion_initialStrain = & - (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**1 / 1. * & - lattice_thermalExpansion33(1:3,1:3,1,phase) + & ! constant coefficient - (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**2 / 2. * & - lattice_thermalExpansion33(1:3,1:3,2,phase) + & ! linear coefficient - (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**3 / 3. * & - lattice_thermalExpansion33(1:3,1:3,3,phase) ! quadratic coefficient + + kinematics_thermal_expansion_initialStrain = & + (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**1 / 1. * & + lattice_thermalExpansion33(1:3,1:3,1,phase) + & ! constant coefficient + (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**2 / 2. * & + lattice_thermalExpansion33(1:3,1:3,2,phase) + & ! linear coefficient + (temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**3 / 3. * & + lattice_thermalExpansion33(1:3,1:3,3,phase) ! quadratic coefficient end function kinematics_thermal_expansion_initialStrain + !-------------------------------------------------------------------------------------------------- !> @brief contains the constitutive equation for calculating the velocity gradient !-------------------------------------------------------------------------------------------------- subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar, ipc, ip, el) - use material, only: & - material_phase, & - material_homog, & - temperature, & - temperatureRate, & - thermalMapping - use lattice, only: & - lattice_thermalExpansion33, & - lattice_referenceTemperature - - implicit none - integer(pInt), intent(in) :: & - ipc, & !< grain number - ip, & !< integration point number - el !< element number - real(pReal), intent(out), dimension(3,3) :: & - Li !< thermal velocity gradient - real(pReal), intent(out), dimension(3,3,3,3) :: & - dLi_dTstar !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero) - integer(pInt) :: & - phase, & - homog, offset - real(pReal) :: & - T, TRef, TDot - - phase = material_phase(ipc,ip,el) - homog = material_homog(ip,el) - offset = thermalMapping(homog)%p(ip,el) - T = temperature(homog)%p(offset) - TDot = temperatureRate(homog)%p(offset) - TRef = lattice_referenceTemperature(phase) - - Li = TDot * ( & - lattice_thermalExpansion33(1:3,1:3,1,phase)*(T - TRef)**0 & ! constant coefficient - + lattice_thermalExpansion33(1:3,1:3,2,phase)*(T - TRef)**1 & ! linear coefficient - + lattice_thermalExpansion33(1:3,1:3,3,phase)*(T - TRef)**2 & ! quadratic coefficient - ) / & - (1.0_pReal & - + lattice_thermalExpansion33(1:3,1:3,1,phase)*(T - TRef)**1 / 1. & - + lattice_thermalExpansion33(1:3,1:3,2,phase)*(T - TRef)**2 / 2. & - + lattice_thermalExpansion33(1:3,1:3,3,phase)*(T - TRef)**3 / 3. & - ) - dLi_dTstar = 0.0_pReal + use material, only: & + material_phase, & + material_homog, & + temperature, & + temperatureRate, & + thermalMapping + use lattice, only: & + lattice_thermalExpansion33, & + lattice_referenceTemperature + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), intent(out), dimension(3,3) :: & + Li !< thermal velocity gradient + real(pReal), intent(out), dimension(3,3,3,3) :: & + dLi_dTstar !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero) + integer(pInt) :: & + phase, & + homog, offset + real(pReal) :: & + T, TRef, TDot + + phase = material_phase(ipc,ip,el) + homog = material_homog(ip,el) + offset = thermalMapping(homog)%p(ip,el) + T = temperature(homog)%p(offset) + TDot = temperatureRate(homog)%p(offset) + TRef = lattice_referenceTemperature(phase) + + Li = TDot * ( & + lattice_thermalExpansion33(1:3,1:3,1,phase)*(T - TRef)**0 & ! constant coefficient + + lattice_thermalExpansion33(1:3,1:3,2,phase)*(T - TRef)**1 & ! linear coefficient + + lattice_thermalExpansion33(1:3,1:3,3,phase)*(T - TRef)**2 & ! quadratic coefficient + ) / & + (1.0_pReal & + + lattice_thermalExpansion33(1:3,1:3,1,phase)*(T - TRef)**1 / 1. & + + lattice_thermalExpansion33(1:3,1:3,2,phase)*(T - TRef)**2 / 2. & + + lattice_thermalExpansion33(1:3,1:3,3,phase)*(T - TRef)**3 / 3. & + ) + dLi_dTstar = 0.0_pReal end subroutine kinematics_thermal_expansion_LiAndItsTangent diff --git a/src/material.f90 b/src/material.f90 index 8a8f36a55..291b73910 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -921,7 +921,7 @@ subroutine material_allocatePlasticState(phase,NofMyPhase,& sizeState,sizeDotState,sizeDeltaState,& Nslip,Ntwin,Ntrans) use numerics, only: & - numerics_integrator2 => numerics_integrator ! compatibility hack + numerics_integrator implicit none integer(pInt), intent(in) :: & @@ -933,8 +933,6 @@ subroutine material_allocatePlasticState(phase,NofMyPhase,& Nslip, & Ntwin, & Ntrans - integer(pInt) :: numerics_integrator ! compatibility hack - numerics_integrator = numerics_integrator2(1) ! compatibility hack plasticState(phase)%sizeState = sizeState plasticState(phase)%sizeDotState = sizeDotState @@ -971,7 +969,7 @@ end subroutine material_allocatePlasticState subroutine material_allocateSourceState(phase,of,NofMyPhase,& sizeState,sizeDotState,sizeDeltaState) use numerics, only: & - numerics_integrator2 => numerics_integrator ! compatibility hack + numerics_integrator implicit none integer(pInt), intent(in) :: & @@ -979,8 +977,6 @@ subroutine material_allocateSourceState(phase,of,NofMyPhase,& of, & NofMyPhase, & sizeState, sizeDotState,sizeDeltaState - integer(pInt) :: numerics_integrator ! compatibility hack - numerics_integrator = numerics_integrator2(1) ! compatibility hack sourceState(phase)%p(of)%sizeState = sizeState sourceState(phase)%p(of)%sizeDotState = sizeDotState diff --git a/src/numerics.f90 b/src/numerics.f90 index 1678d0c48..9727a04a7 100644 --- a/src/numerics.f90 +++ b/src/numerics.f90 @@ -23,12 +23,10 @@ module numerics pert_method = 1_pInt, & !< method used in perturbation technique for tangent randomSeed = 0_pInt, & !< fixed seeding for pseudo-random number generator, Default 0: use random seed worldrank = 0_pInt, & !< MPI worldrank (/=0 for MPI simulations only) - worldsize = 0_pInt !< MPI worldsize (/=0 for MPI simulations only) + worldsize = 0_pInt, & !< MPI worldsize (/=0 for MPI simulations only) + numerics_integrator = 1_pInt !< method used for state integration Default 1: fix-point iteration integer(4), protected, public :: & DAMASK_NumThreadsInt = 0 !< value stored in environment variable DAMASK_NUM_THREADS, set to zero if no OpenMP directive - !< ToDo: numerics_integrator is an array for historical reasons, only element 1 is used! - integer(pInt), dimension(2), protected, public :: & - numerics_integrator = 1_pInt !< method used for state integration (central & perturbed state), Default 1: fix-point iteration for both states real(pReal), protected, public :: & relevantStrain = 1.0e-7_pReal, & !< strain increment considered significant (used by crystallite to determine whether strain inc is considered significant) defgradTolerance = 1.0e-7_pReal, & !< deviation of deformation gradient that is still allowed (used by CPFEM to determine outdated ffn1) @@ -466,7 +464,7 @@ subroutine numerics_init write(6,'(a24,1x,es8.1)') ' rTol_crystalliteState: ',rTol_crystalliteState write(6,'(a24,1x,es8.1)') ' rTol_crystalliteStress: ',rTol_crystalliteStress write(6,'(a24,1x,es8.1)') ' aTol_crystalliteStress: ',aTol_crystalliteStress - write(6,'(a24,2(1x,i8))') ' integrator: ',numerics_integrator + write(6,'(a24,1x,i8)') ' integrator: ',numerics_integrator write(6,'(a24,1x,L8)') ' use ping pong scheme: ',usepingpong write(6,'(a24,1x,es8.1,/)')' unitlength: ',numerics_unitlength @@ -589,7 +587,7 @@ subroutine numerics_init if (rTol_crystalliteState <= 0.0_pReal) call IO_error(301_pInt,ext_msg='rTol_crystalliteState') if (rTol_crystalliteStress <= 0.0_pReal) call IO_error(301_pInt,ext_msg='rTol_crystalliteStress') if (aTol_crystalliteStress <= 0.0_pReal) call IO_error(301_pInt,ext_msg='aTol_crystalliteStress') - if (any(numerics_integrator <= 0_pInt) .or. any(numerics_integrator >= 6_pInt)) & + if (numerics_integrator <= 0_pInt .or. numerics_integrator >= 6_pInt) & call IO_error(301_pInt,ext_msg='integrator') if (numerics_unitlength <= 0.0_pReal) call IO_error(301_pInt,ext_msg='unitlength') if (absTol_RGC <= 0.0_pReal) call IO_error(301_pInt,ext_msg='absTol_RGC') diff --git a/src/source_thermal_dissipation.f90 b/src/source_thermal_dissipation.f90 index ef0843239..026a71726 100644 --- a/src/source_thermal_dissipation.f90 +++ b/src/source_thermal_dissipation.f90 @@ -1,4 +1,5 @@ !-------------------------------------------------------------------------------------------------- +!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH !> @brief material subroutine for thermal source due to plastic dissipation !> @details to be done diff --git a/src/source_thermal_externalheat.f90 b/src/source_thermal_externalheat.f90 index 724d36f17..2bf4cac9c 100644 --- a/src/source_thermal_externalheat.f90 +++ b/src/source_thermal_externalheat.f90 @@ -1,53 +1,44 @@ !-------------------------------------------------------------------------------------------------- +!> @author Martin Diehl, Max-Planck-Institut für Eisenforschung GmbH !> @author Pratheek Shanthraj, Max-Planck-Institut für Eisenforschung GmbH !> @author Philip Eisenlohr, Michigan State University !> @brief material subroutine for variable heat source -!> @details to be done !-------------------------------------------------------------------------------------------------- module source_thermal_externalheat - use prec, only: & - pReal, & - pInt + use prec, only: & + pReal, & + pInt - implicit none - private - integer(pInt), dimension(:), allocatable, public, protected :: & - source_thermal_externalheat_sizePostResults, & !< cumulative size of post results - source_thermal_externalheat_offset, & !< which source is my current thermal dissipation mechanism? - source_thermal_externalheat_instance !< instance of thermal dissipation source mechanism + implicit none + private + integer(pInt), dimension(:), allocatable, public, protected :: & + source_thermal_externalheat_offset, & !< which source is my current thermal dissipation mechanism? + source_thermal_externalheat_instance !< instance of thermal dissipation source mechanism - integer(pInt), dimension(:,:), allocatable, target, public :: & - source_thermal_externalheat_sizePostResult !< size of each post result output + integer(pInt), dimension(:,:), allocatable, target, public :: & + source_thermal_externalheat_sizePostResult !< size of each post result output - character(len=64), dimension(:,:), allocatable, target, public :: & - source_thermal_externalheat_output !< name of each post result output + character(len=64), dimension(:,:), allocatable, target, public :: & + source_thermal_externalheat_output !< name of each post result output - integer(pInt), dimension(:), allocatable, target, public :: & - source_thermal_externalheat_Noutput !< number of outputs per instance of this source + integer(pInt), dimension(:), allocatable, target, public :: & + source_thermal_externalheat_Noutput !< number of outputs per instance of this source - integer(pInt), dimension(:), allocatable, private :: & - source_thermal_externalheat_nIntervals + type, private :: tParameters !< container type for internal constitutive parameters + real(pReal), dimension(:), allocatable :: & + time, & + heat_rate + integer(pInt) :: & + nIntervals + end type tParameters - real(pReal), dimension(:,:), allocatable, private :: & - source_thermal_externalheat_time, & - source_thermal_externalheat_rate + type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) - type, private :: tParameters !< container type for internal constitutive parameters - real(pReal), dimension(:), allocatable :: & - time, & - rate - integer(pInt) :: & - nInterval - end type tParameters - - type(tParameters), dimension(:), allocatable, private :: param !< containers of constitutive parameters (len Ninstance) - - - public :: & - source_thermal_externalheat_init, & - source_thermal_externalheat_dotState, & - source_thermal_externalheat_getRateAndItsTangent + public :: & + source_thermal_externalheat_init, & + source_thermal_externalheat_dotState, & + source_thermal_externalheat_getRateAndItsTangent contains @@ -56,143 +47,74 @@ contains !> @brief module initialization !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- -subroutine source_thermal_externalheat_init(fileUnit) - use debug, only: & - debug_level,& - debug_constitutive,& - debug_levelBasic - use IO, only: & - IO_read, & - IO_lc, & - IO_getTag, & - IO_isBlank, & - IO_stringPos, & - IO_stringValue, & - IO_floatValue, & - IO_intValue, & - IO_warning, & - IO_error, & - IO_EOF - use material, only: & - material_allocateSourceState, & - phase_source, & - phase_Nsources, & - phase_Noutput, & - SOURCE_thermal_externalheat_label, & - SOURCE_thermal_externalheat_ID, & - material_phase, & - sourceState - use config, only: & - config_phase, & - material_Nphase, & - MATERIAL_partPhase - - implicit none - integer(pInt), intent(in) :: fileUnit - - integer(pInt), allocatable, dimension(:) :: chunkPos - integer(pInt) :: maxNinstance,phase,instance,source,sourceOffset - integer(pInt) :: sizeState, sizeDotState, sizeDeltaState - integer(pInt) :: NofMyPhase,interval,p - character(len=65536) :: & - tag = '', & - line = '' - real(pReal), allocatable, dimension(:,:) :: temp_time, temp_rate - - write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_externalheat_label//' init -+>>>' - +subroutine source_thermal_externalheat_init + use debug, only: & + debug_level,& + debug_constitutive,& + debug_levelBasic + use material, only: & + material_allocateSourceState, & + material_phase, & + phase_source, & + phase_Nsources, & + phase_Noutput, & + SOURCE_thermal_externalheat_label, & + SOURCE_thermal_externalheat_ID + use config, only: & + config_phase, & + material_Nphase, & + MATERIAL_partPhase - maxNinstance = int(count(phase_source == SOURCE_thermal_externalheat_ID),pInt) - if (maxNinstance == 0_pInt) return - if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & - write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance + implicit none - allocate(source_thermal_externalheat_offset(material_Nphase), source=0_pInt) - allocate(source_thermal_externalheat_instance(material_Nphase), source=0_pInt) - do phase = 1, material_Nphase - source_thermal_externalheat_instance(phase) = count(phase_source(:,1:phase) == SOURCE_thermal_externalheat_ID) - do source = 1, phase_Nsources(phase) - if (phase_source(source,phase) == SOURCE_thermal_externalheat_ID) & - source_thermal_externalheat_offset(phase) = source - enddo - enddo - - allocate(source_thermal_externalheat_sizePostResults(maxNinstance), source=0_pInt) - allocate(source_thermal_externalheat_sizePostResult(maxval(phase_Noutput),maxNinstance),source=0_pInt) - allocate(source_thermal_externalheat_output (maxval(phase_Noutput),maxNinstance)) - source_thermal_externalheat_output = '' - allocate(source_thermal_externalheat_Noutput(maxNinstance), source=0_pInt) - - allocate(source_thermal_externalheat_nIntervals(maxNinstance), source=0_pInt) - allocate(temp_time(maxNinstance,1000), source=0.0_pReal) - allocate(temp_rate(maxNinstance,1000), source=0.0_pReal) - - do p=1, size(config_phase) - if (all(phase_source(:,p) /= SOURCE_thermal_externalheat_ID)) cycle - enddo - - rewind(fileUnit) - phase = 0_pInt - do while (trim(line) /= IO_EOF .and. IO_lc(IO_getTag(line,'<','>')) /= MATERIAL_partPhase) ! wind forward to - line = IO_read(fileUnit) - enddo + real(pReal), allocatable, dimension(:) :: tempVar + integer(pInt) :: maxNinstance,instance,source,sourceOffset + integer(pInt) :: NofMyPhase,p - parsingFile: do while (trim(line) /= IO_EOF) ! read through sections of phase part - line = IO_read(fileUnit) - if (IO_isBlank(line)) cycle ! skip empty lines - if (IO_getTag(line,'<','>') /= '') then ! stop at next part - line = IO_read(fileUnit, .true.) ! reset IO_read - exit - endif - if (IO_getTag(line,'[',']') /= '') then ! next phase section - phase = phase + 1_pInt ! advance phase section counter - cycle ! skip to next line - endif - - if (phase > 0_pInt ) then; if (any(phase_source(:,phase) == SOURCE_thermal_externalheat_ID)) then ! do not short-circuit here (.and. with next if statemen). It's not safe in Fortran - - instance = source_thermal_externalheat_instance(phase) ! which instance of my source is present phase - chunkPos = IO_stringPos(line) - tag = IO_lc(IO_stringValue(line,chunkPos,1_pInt)) ! extract key - select case(tag) - case ('externalheat_time','externalheat_rate') - if (chunkPos(1) <= 2_pInt) & - call IO_error(150_pInt,ext_msg=trim(tag)//' ('//SOURCE_thermal_externalheat_label//')') - if ( source_thermal_externalheat_nIntervals(instance) > 0_pInt & - .and. source_thermal_externalheat_nIntervals(instance) /= chunkPos(1) - 2_pInt) & - call IO_error(150_pInt,ext_msg=trim(tag)//' ('//SOURCE_thermal_externalheat_label//')') - - source_thermal_externalheat_nIntervals(instance) = chunkPos(1) - 2_pInt - do interval = 1, source_thermal_externalheat_nIntervals(instance) + 1_pInt - select case(tag) - case ('externalheat_time') - temp_time(instance, interval) = IO_floatValue(line,chunkPos,1_pInt + interval) - case ('externalheat_rate') - temp_rate(instance, interval) = IO_floatValue(line,chunkPos,1_pInt + interval) - end select - enddo - end select - endif; endif - enddo parsingFile + write(6,'(/,a)') ' <<<+- source_'//SOURCE_thermal_externalheat_label//' init -+>>>' - allocate(source_thermal_externalheat_time(maxNinstance,maxval(source_thermal_externalheat_nIntervals)+1_pInt), source=0.0_pReal) - allocate(source_thermal_externalheat_rate(maxNinstance,maxval(source_thermal_externalheat_nIntervals)+1_pInt), source=0.0_pReal) - - initializeInstances: do phase = 1_pInt, material_Nphase - if (any(phase_source(:,phase) == SOURCE_thermal_externalheat_ID)) then - NofMyPhase = count(material_phase==phase) - instance = source_thermal_externalheat_instance(phase) - sourceOffset = source_thermal_externalheat_offset(phase) - source_thermal_externalheat_time(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) = & - temp_time(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) - source_thermal_externalheat_rate(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) = & - temp_rate(instance,1:source_thermal_externalheat_nIntervals(instance)+1_pInt) - - call material_allocateSourceState(phase,sourceOffset,NofMyPhase,1_pInt,1_pInt,0_pInt) + + maxNinstance = int(count(phase_source == SOURCE_thermal_externalheat_ID),pInt) + if (maxNinstance == 0_pInt) return + if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0_pInt) & + write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance + + allocate(source_thermal_externalheat_offset(material_Nphase), source=0_pInt) + allocate(source_thermal_externalheat_instance(material_Nphase), source=0_pInt) + + do p = 1, material_Nphase + source_thermal_externalheat_instance(p) = count(phase_source(:,1:p) == SOURCE_thermal_externalheat_ID) + do source = 1, phase_Nsources(p) + if (phase_source(source,p) == SOURCE_thermal_externalheat_ID) & + source_thermal_externalheat_offset(p) = source + enddo + enddo - endif + allocate(source_thermal_externalheat_sizePostResult(maxval(phase_Noutput),maxNinstance),source=0_pInt) + allocate(source_thermal_externalheat_output (maxval(phase_Noutput),maxNinstance)) + source_thermal_externalheat_output = '' + allocate(source_thermal_externalheat_Noutput(maxNinstance), source=0_pInt) - enddo initializeInstances + allocate(param(maxNinstance)) + + do p=1, size(config_phase) + if (all(phase_source(:,p) /= SOURCE_thermal_externalheat_ID)) cycle + instance = source_thermal_externalheat_instance(p) + sourceOffset = source_thermal_externalheat_offset(p) + NofMyPhase=count(material_phase==p) + + tempVar = config_phase(p)%getFloats('externalheat_time') + param(instance)%nIntervals = size(tempVar) - 1_pInt + + param(instance)%time= tempVar + + tempVar = config_phase(p)%getFloats('externalheat_rate',requiredSize = size(tempVar)) + param(instance)%heat_rate = tempVar + + call material_allocateSourceState(p,sourceOffset,NofMyPhase,1_pInt,1_pInt,0_pInt) + + enddo + end subroutine source_thermal_externalheat_init @@ -245,16 +167,16 @@ subroutine source_thermal_externalheat_getRateAndItsTangent(TDot, dTDot_dT, phas instance = source_thermal_externalheat_instance(phase) sourceOffset = source_thermal_externalheat_offset(phase) - do interval = 1, source_thermal_externalheat_nIntervals(instance) ! scan through all rate segments + do interval = 1, param(instance)%nIntervals ! scan through all rate segments frac_time = (sourceState(phase)%p(sourceOffset)%state(1,constituent) - & - source_thermal_externalheat_time(instance,interval)) / & - (source_thermal_externalheat_time(instance,interval+1) - & - source_thermal_externalheat_time(instance,interval)) ! fractional time within segment + param(instance)%time(interval)) / & + (param(instance)%time(interval+1) - & + param(instance)%time(interval)) ! fractional time within segment if ( (frac_time < 0.0_pReal .and. interval == 1) & - .or. (frac_time >= 1.0_pReal .and. interval == source_thermal_externalheat_nIntervals(instance)) & + .or. (frac_time >= 1.0_pReal .and. interval == param(instance)%nIntervals) & .or. (frac_time >= 0.0_pReal .and. frac_time < 1.0_pReal) ) & - TDot = source_thermal_externalheat_rate(instance,interval ) * (1.0_pReal - frac_time) + & - source_thermal_externalheat_rate(instance,interval+1) * frac_time ! interpolate heat rate between segment boundaries... + TDot = param(instance)%heat_rate(interval ) * (1.0_pReal - frac_time) + & + param(instance)%heat_rate(interval+1) * frac_time ! interpolate heat rate between segment boundaries... ! ...or extrapolate if outside of bounds enddo dTDot_dT = 0.0 diff --git a/src/thermal_adiabatic.f90 b/src/thermal_adiabatic.f90 index 947e28777..c3290bdfe 100644 --- a/src/thermal_adiabatic.f90 +++ b/src/thermal_adiabatic.f90 @@ -43,8 +43,6 @@ contains !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- subroutine thermal_adiabatic_init - use IO, only: & - IO_error use material, only: & thermal_type, & thermal_typeInstance, & diff --git a/src/thermal_conduction.f90 b/src/thermal_conduction.f90 index 3d5ca892e..88da0529b 100644 --- a/src/thermal_conduction.f90 +++ b/src/thermal_conduction.f90 @@ -44,8 +44,6 @@ contains !> @details reads in material parameters, allocates arrays, and does sanity checks !-------------------------------------------------------------------------------------------------- subroutine thermal_conduction_init - use IO, only: & - IO_error use material, only: & thermal_type, & thermal_typeInstance, &