Hi,
I'm currently working on the modelling of a thermodynamic process.
Briefly, I cool down a solution (water + polymer) from -5°C to -15°C to induce a phase separation. At the end (and after removing of the water by lyophilisation) I obtain a porous sponge like material.
The process uses a home made cooling system which can be described like this:
- A Peltier module
- An aluminium layer recovered by teflon (And also a layer of ethanol)
- The mold (in stainless steel) in which I pour the (polymer+water) solution
- A layer of air.
The bottom part of the Peltier is kept at -15°C during the all experiment. Then the top part Peltier (x=0) heat to -5°C the full system until thermal equilibrium. At t0, the Peltier lets the temperature in x=0 goes down linearly to -15°C in a controlled way by the adjustement of the cooling rate. A thermocouple measure the temperature in the aluminium part and at ~-10°C (in the aluminium) I observe an increase of the measured temperature due to the phase change of the solution.
I try to model the full process to determine how long is necessary to reach -15°C in the (polymer+water) solution. At the moment I already made it but without the introduction of the phase change of the solution which can be considered as a heat source. Normally I should obtain a decreasing surface with a kind of 2D gaussian in x= polymer solution... But It doesn't work...
At the moment, I use this method:
Definition of the variables
Cooling rate of the Peltier module (in °C.min) : R1 := .5; (As T(0,t)=T0-R1*time)
Cooling rate of the air in the cooling chamber (in °C.min) : R2 := .95*R1; Thid should represent the convection function which is not represented here (It's obviously an approximation) as I don't know how to work with real convection... I need it as boundary conditions.
Initial time of cooling (in min): t0 := 20; (After 20 min I switch on the system and it begins to cool down)
Initial temperature of the system (in °C) : T0 := -5;
Final temperature of the system (in °C) : Tf := -15;
Cooling time in x=0 (in min): tau1 := (T0-Tf)/R1;
Cooling time of air, x= top of the full system (in min): tau2 := (T0-Tf)/R2;
Time of experiment (in min): texp := t0+tau1+50; I wait ~ 50minutes at the end, to reach the equilibrium but this value is completely arbitrariness.
Description of the geometry
1. Thickness (in m):
Aluminium x[1] := 0.6e-2;
Teflon film x[2] := 0.1e-3;
Ethanol x[3] := 0.1e-3;
Stainless steel x[4] := 0.2e-2;
Polymer solution x[5] := 0.2e-2;
Air x[6] := .1;
Size of the system (in m): Ltot := sum(x[n], n = 1 .. 6);
Region of interest (in m): Lint := sum(x[n], n = 1 .. 5);
2. Values of the thermal diffusivity in function of x
Thermal diffusivity (in m2.s-1)
Aluminium alpha1 := 1/(8.418*10^(-5));
Teflon film alpha2 := 1/(.124*10^(-6));
Ethanol alpha3 := 1/(7*10^(-8));
Stainless steel alpha4 := 1/(3.428*10^(-6));
polymer solution (Assimilated to water) alpha5 := 1/(.143*10^(-6));
Air alpha6 := 1/(1.9*10^(-5));
Functions for boundary conditions
1. Thermal diffusivity function (Considered as a constant on this temperature range)
alpha := (x) -> piecewise(0 < x and x
2. Boundary conditions in x = 0
Ti := (t)-> piecewise(0 < t and t
3. Boundary conditions in x = Ltot (The air cools down slowly but we choose to introduce it as a boundary condition)
Tl := (t) -> piecewise(0 < t and t
Solving of the heat equation
1. Initial and boundary conditions
IBC := {T(0, t) = Ti(t), T(x, 0) = T0, T(Ltot, t) = Tl(t)};
2. Numeric solving of the heat equation
PDE := diff(T(x, t), x, x)-alpha(x)*(diff(T(x, t), t)) = 0;
pds := pdsolve(PDE, IBC, numeric);
3. Plotting domain
P1 := pds:-plot3d(x = 0 .. Lint, t = 0 .. texp);
4. Plot
plots[display]({P1}, axes = boxed);
This gives a nice estimation of the temperature evolution in the system, I can change a lot of parameters for optimisation and I see directly the temperature evolution in the system according to the cooling rate in x=0.
HOWEVER, I don't know how to introduce the local temperature increase during the phase change, I tried to add a third boundary condition to solve it by introducing a 2D gaussian curve but the software refuses this and it's not really the best I think. The best would be to add a source term in the heat equation but I don't know how... SOMEBODY TO HELP ME ?
I don't know if I'm clear enough but it would be really nice to help me on this topic which is, of course, fascinating... No ?
Thanks in advance,
Lisa
