From 55b88e47b711e189c729aeb3b9ec45aed3300d26 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 9 Jan 2013 18:08:08 +0000 Subject: [PATCH] fixed wrong temperature when using spectral solver --- code/DAMASK_spectral_driver.f90 | 32 +++++++++++------------ code/DAMASK_spectral_solverAL.f90 | 12 +++++---- code/DAMASK_spectral_solverBasic.f90 | 17 ++++++------ code/DAMASK_spectral_solverBasicPETSc.f90 | 9 ++++--- 4 files changed, 38 insertions(+), 32 deletions(-) diff --git a/code/DAMASK_spectral_driver.f90 b/code/DAMASK_spectral_driver.f90 index ee3c1c35e..af0658e85 100644 --- a/code/DAMASK_spectral_driver.f90 +++ b/code/DAMASK_spectral_driver.f90 @@ -293,12 +293,12 @@ program DAMASK_spectral_Driver ! doing initialization depending on selected solver select case (myspectralsolver) case (DAMASK_spectral_SolverBasic_label) - call basic_init() + call basic_init(loadCases(1)%temperature) #ifdef PETSc case (DAMASK_spectral_SolverBasicPETSc_label) - call basicPETSc_init() + call basicPETSc_init(loadCases(1)%temperature) case (DAMASK_spectral_SolverAL_label) - call AL_init() + call AL_init(loadCases(1)%temperature) #endif case default call IO_error(error_ID = 891, ext_msg = trim(myspectralsolver)) @@ -398,27 +398,27 @@ program DAMASK_spectral_Driver ! calculate solution case (DAMASK_spectral_SolverBasic_label) solres = basic_solution (& - incInfo, guess,timeinc,timeinc_old, & - P_BC = loadCases(currentLoadCase)%P, & - F_BC = loadCases(currentLoadCase)%deformation, & - temperature_bc = loadCases(currentLoadCase)%temperature, & - rotation_BC = loadCases(currentLoadCase)%rotation) + incInfo, guess,timeinc,timeinc_old, & + P_BC = loadCases(currentLoadCase)%P, & + F_BC = loadCases(currentLoadCase)%deformation, & + temperature_bc = loadCases(currentLoadCase)%temperature, & + rotation_BC = loadCases(currentLoadCase)%rotation) #ifdef PETSc case (DAMASK_spectral_SolverBasicPETSC_label) solres = BasicPETSC_solution (& incInfo, guess,timeinc,timeinc_old, & - P_BC = loadCases(currentLoadCase)%P, & - F_BC = loadCases(currentLoadCase)%deformation, & - temperature_bc = loadCases(currentLoadCase)%temperature, & - rotation_BC = loadCases(currentLoadCase)%rotation) + P_BC = loadCases(currentLoadCase)%P, & + F_BC = loadCases(currentLoadCase)%deformation, & + temperature_bc = loadCases(currentLoadCase)%temperature, & + rotation_BC = loadCases(currentLoadCase)%rotation) case (DAMASK_spectral_SolverAL_label) solres = AL_solution (& incInfo, guess,timeinc,timeinc_old, & - P_BC = loadCases(currentLoadCase)%P, & - F_BC = loadCases(currentLoadCase)%deformation, & - temperature_bc = loadCases(currentLoadCase)%temperature, & - rotation_BC = loadCases(currentLoadCase)%rotation) + P_BC = loadCases(currentLoadCase)%P, & + F_BC = loadCases(currentLoadCase)%deformation, & + temperature_bc = loadCases(currentLoadCase)%temperature, & + rotation_BC = loadCases(currentLoadCase)%rotation) #endif end select !-------------------------------------------------------------------------------------------------- diff --git a/code/DAMASK_spectral_solverAL.f90 b/code/DAMASK_spectral_solverAL.f90 index 2464c4198..d62caab17 100644 --- a/code/DAMASK_spectral_solverAL.f90 +++ b/code/DAMASK_spectral_solverAL.f90 @@ -28,6 +28,7 @@ module DAMASK_spectral_solverAL type tSolutionParams real(pReal), dimension(3,3) :: P_BC, rotation_BC real(pReal) :: timeinc + real(pReal) :: temperature end type tSolutionParams type(tSolutionParams), private :: params @@ -46,9 +47,7 @@ module DAMASK_spectral_solverAL F_lambda_lastInc, & !< field of previous incompatible deformation gradient Fdot, & !< field of assumed rate of compatible deformation gradient F_lambdaDot !< field of assumed rate of incopatible deformation gradient - real(pReal), private :: & - temperature !< temperature, no spatial quantity at the moment - + !-------------------------------------------------------------------------------------------------- ! stress, stiffness and compliance average etc. real(pReal), private, dimension(3,3) :: & @@ -75,7 +74,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields and fills them with data, potentially from restart info !-------------------------------------------------------------------------------------------------- -subroutine AL_init() +subroutine AL_init(temperature) use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran >4.6 at the moment) use IO, only: & IO_read_JobBinaryFile, & @@ -101,6 +100,8 @@ subroutine AL_init() math_invSym3333 implicit none + real(pReal), intent(inout) :: & + temperature #include #include integer(pInt) :: i,j,k @@ -344,6 +345,7 @@ type(tSolutionState) function & params%P_BC = P_BC%values params%rotation_BC = rotation_BC params%timeinc = timeinc + params%temperature = temperature_bc call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr) CHKERRQ(ierr) @@ -429,7 +431,7 @@ subroutine AL_formResidual(in,x_scal,f_scal,dummy,ierr) !-------------------------------------------------------------------------------------------------- ! evaluate constitutive response - call Utilities_constitutiveResponse(F_lastInc,F,temperature,params%timeinc, & + call Utilities_constitutiveResponse(F_lastInc,F,params%temperature,params%timeinc, & residual_F,C,P_av,ForwardData,params%rotation_BC) ForwardData = .False. diff --git a/code/DAMASK_spectral_solverBasic.f90 b/code/DAMASK_spectral_solverBasic.f90 index a8768a921..9be1c1808 100644 --- a/code/DAMASK_spectral_solverBasic.f90 +++ b/code/DAMASK_spectral_solverBasic.f90 @@ -27,8 +27,7 @@ module DAMASK_spectral_SolverBasic F, & !< deformation gradient field F_lastInc, & !< deformation gradient field last increment Fdot !< assumed rate for F n to F n+1 - real(pReal), private :: temperature !< not pointwise - + !-------------------------------------------------------------------------------------------------- ! stress, stiffness and compliance average etc. real(pReal), private, dimension(3,3) :: & @@ -44,7 +43,7 @@ contains !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields and fills them with data, potentially from restart info !-------------------------------------------------------------------------------------------------- -subroutine basic_init() +subroutine basic_init(temperature) use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran >4.6 at the moment) use IO, only: & @@ -69,6 +68,8 @@ subroutine basic_init() mesh_deformedCoordsFFT implicit none + real(pReal), intent(inout) :: & + temperature real(pReal), dimension(3,3,res(1),res(2),res(3)) :: P integer(pInt) :: & i, j, k @@ -124,7 +125,7 @@ subroutine basic_init() endif mesh_ipCoordinates = 0.0_pReal !reshape(mesh_deformedCoordsFFT(geomdim,& !reshape(F,[3,3,res(1),res(2),res(3)])),[3,1,mesh_NcpElems]) - call Utilities_constitutiveResponse(F,F,temperature,0.0_pReal,P,C,temp33_Real,.false.,math_I3) ! constitutive response with no deformation in no time to get reference stiffness + call Utilities_constitutiveResponse(F,F,temperature,0.0_pReal,P,C,temp33_Real,.false.,math_I3) ! constitutive response with no deformation in no time to get reference stiffness if (restartInc == 1_pInt) then ! use initial stiffness as reference stiffness temp3333_Real = C endif @@ -183,8 +184,7 @@ type(tSolutionState) function & ! input data for solution real(pReal), intent(in) :: & timeinc, & !< increment in time for current solution - timeinc_old, & !< increment in time of last increment - temperature_bc !< temperature (I wonder, why it is intent(in) + timeinc_old !< increment in time of last increment logical, intent(in) :: & guess !< if .false., assume homogeneous addon type(tBoundaryCondition), intent(in) :: & @@ -194,7 +194,8 @@ type(tSolutionState) function & incInfo !< string with information of current increment for output to screen real(pReal), dimension(3,3), intent(in) :: & rotation_BC !< rotation load to lab - + real(pReal), intent(inout) :: & + temperature_bc real(pReal), dimension(3,3,3,3) :: & S !< current average compliance real(pReal), dimension(3,3) :: & @@ -287,7 +288,7 @@ type(tSolutionState) function & ! evaluate constitutive response F_aim_lastIter = F_aim basic_solution%termIll = .false. - call Utilities_constitutiveResponse(F_lastInc,F,temperature,timeinc,& + call Utilities_constitutiveResponse(F_lastInc,F,temperature_bc,timeinc,& P,C,P_av,ForwardData,rotation_BC) basic_solution%termIll = terminallyIll terminallyIll = .false. diff --git a/code/DAMASK_spectral_solverBasicPETSc.f90 b/code/DAMASK_spectral_solverBasicPETSc.f90 index d90b1e8e7..0b3ea970c 100644 --- a/code/DAMASK_spectral_solverBasicPETSc.f90 +++ b/code/DAMASK_spectral_solverBasicPETSc.f90 @@ -28,6 +28,7 @@ module DAMASK_spectral_SolverBasicPETSc type tSolutionParams real(pReal), dimension(3,3) :: P_BC, rotation_BC real(pReal) :: timeinc + real(pReal) :: temperature end type tSolutionParams type(tSolutionParams), private :: params @@ -41,7 +42,6 @@ module DAMASK_spectral_SolverBasicPETSc !-------------------------------------------------------------------------------------------------- ! common pointwise data real(pReal), private, dimension(:,:,:,:,:), allocatable :: F_lastInc, Fdot - real(pReal) :: temperature !-------------------------------------------------------------------------------------------------- ! stress, stiffness and compliance average etc. @@ -69,7 +69,7 @@ module DAMASK_spectral_SolverBasicPETSc !-------------------------------------------------------------------------------------------------- !> @brief allocates all neccessary fields and fills them with data, potentially from restart info !-------------------------------------------------------------------------------------------------- -subroutine basicPETSc_init() +subroutine basicPETSc_init(temperature) use, intrinsic :: iso_fortran_env ! to get compiler_version and compiler_options (at least for gfortran >4.6 at the moment) use IO, only: & IO_read_JobBinaryFile, & @@ -96,6 +96,8 @@ subroutine basicPETSc_init() math_invSym3333 implicit none + real(pReal), intent(inout) :: & + temperature #include #include integer(pInt) :: i,j,k @@ -307,6 +309,7 @@ type(tSolutionState) function & params%P_BC = P_BC%values params%rotation_BC = rotation_BC params%timeinc = timeinc + params%temperature = temperature_BC call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr) CHKERRQ(ierr) @@ -370,7 +373,7 @@ subroutine BasicPETSC_formResidual(myIn,x_scal,f_scal,dummy,ierr) !-------------------------------------------------------------------------------------------------- ! evaluate constitutive response - call Utilities_constitutiveResponse(F_lastInc,x_scal,temperature,params%timeinc, & + call Utilities_constitutiveResponse(F_lastInc,x_scal,params%temperature,params%timeinc, & f_scal,C,P_av,ForwardData,params%rotation_BC) ForwardData = .false.