Discussion Closed This discussion was created more than 6 months ago and has been closed. To start a new discussion with a link back to this one, click here.

Phase transformation with thermal hysteresis

Please login with a confirmed email address before reporting spam

Hi everybody,

my colleagues and I are currently trying to model the austenite-martensite transition in shape memory alloys including the exchange of latent heat. Similar problems (regarding latent heat) have already been discussed here:
www.comsol.com/community/forums/general/thread/31450/
www.comsol.com/community/forums/general/thread/33678
We are using COMSOL 4.3b with the heat transfer and structural mechanics modules, simulating time dependent in 2D or 3D.

We have some problems with the implementation of thermal hysteresis and would be very glad about any help. In the following, I will try to describe the problem as simple as possible (if something is unclear, please ask).

We consider a material with a high temperature phase (austenite) and a low temperature phase (martensite), which I will call A and M in the following. The phase transformation can be triggered either by temperature or by mechanical stress. We use a phenomenologic formula that describes the dependence of the phase fraction Xi of M-phase on temperature and stress (the phase fraction of A-phase is 1-f(T,stress), accordingly):

f(T,stress) = 1/(1+ exp(k0*(T - T0 - c*stress))

This is something like a soft step function, giving values between 0 and 1. The parameters k0, T0 and c are estimated from experimental data; T is the current temperature, for stress we use either solid.mises or solid.tresca. For a start, we have just set the phase fraction Xi = f(T,stress), that means there is one distinct phase fraction at a given temperature and Tresca stress.

We need this phase fraction in several parts of the simulation:
1) The effective Young's modulus is a mixture of the moduli EA and EM of both phases: E = (1-Xi)*EA + Xi*EM
2) During phase transformation, latent heat is absorbed or released. We implemented a local heat source that is proportional to the change of Xi over time:
q = d(Xi, time) * rho * deltaH
(rho: mass density [kg/m³], deltaH: latent heat in [kJ/kg])
You see that the phase fraction depends on stress and temperature and these two, in turn, depend on the phase fraction. Hence we have a fully coupled problem and had to build in some weak constraints to make the simulation run.

So far, we have implemented the model and it works properly for simple geometries. If we start with a block of material in phase A and ramp up the strain in x-direction at constant ambient temperature, the phase fraction Xi will go from 0 to 1 above a certain critical stress and the sample will heat up due to latent heat. (See attached plots stress.png and phasefraction.png)

However, the real material behavior includes hysteresis and this is where our problem lies. There are actually two functions f_AM(T,stress) and f_MA(T,stress), one describing the A->M transformation, and the other the reverse process (M->A), respectively. The parameters k0, T0 and c are not the same in both cases. Also, the material behavior becomes history-dependent. If we again consider a sample of 100% A-phase and ramp up the strain linearly in time, the M phase fraction will increase according to f_AM(T,mises). We can stop (for example) when only 25% of phase A is left (Xi=0.75) and ramp the strain down again. Now, only 75% of the material can take part in the M->A back transformation (described by f_MA(T,mises)), because 25% already are in A-phase. (Compare attached file Graph1.png)

Expressed in pseudo-code, what we want to do is the following:

IF solid.tresca > prev(solid.tresca)
// loading -> Xi may only increase
Xi = 1 - (1-prev(Xi)) * (1-f_AM(T,solid.tresca));
ELSE IF solid.tresca < prev(solid.tresca)
// unloading -> Xi may only decrease
Xi = prev(Xi) * f_MA(T,solid.tresca);
ELSE
// phase fraction remains the same
Xi = prev(Xi);
END

Unfortunately we can not use COMSOL's prev() operator, because we need the solid mechanics module which does not seem to support time-discrete solving. Does anyone have an idea, how to implement such a variable Xi (which depends on the values of itself and other quantities at previous times) for the time-dependent solver?

Thank you very much for reading all this.


Kind regards,
H. Ossmer











1 Reply Last Post 1 oct. 2015, 08:27 UTC−4
COMSOL Moderator

Hello Hinnerk Oßmer

Your Discussion has gone 30 days without a reply. If you still need help with COMSOL and have an on-subscription license, please visit our Support Center for help.

If you do not hold an on-subscription license, you may find an answer in another Discussion or in the Knowledge Base.


Please login with a confirmed email address before reporting spam

Posted: 9 years ago 1 oct. 2015, 08:27 UTC−4
Is there any solution for this Problem yet?
Is there any solution for this Problem yet?

Note that while COMSOL employees may participate in the discussion forum, COMSOL® software users who are on-subscription should submit their questions via the Support Center for a more comprehensive response from the Technical Support team.