0001 function md=icethermal(md,analysis_type)
0002
0003
0004
0005
0006
0007
0008 iceglobal
0009
0010
0011 if strcmpi(md.cluster,'yes'), cluster=1; else cluster=0;end;
0012
0013
0014 if cluster,
0015 error('icethermal error message: parallel support not implemented yet');
0016 end
0017
0018
0019 m_t=CreateFemModel(md,analysis_type);
0020 if strcmpi(analysis_type,'thermalsteady'),
0021 m_m=CreateFemModel(md,'meltingsteady');
0022 else
0023 m_m=CreateFemModel(md,'meltingtransient');
0024 end
0025
0026 velocity=zeros(m_t.gridset.gsize,1);
0027 velocity(1:6:m_t.gridset.gsize)=md.vx(m_t.part)/md.yts;
0028 velocity(2:6:m_t.gridset.gsize)=md.vy(m_t.part)/md.yts;
0029 velocity(3:6:m_t.gridset.gsize)=md.vz(m_t.part)/md.yts;
0030
0031 pressure=zeros(m_t.gridset.gsize,1);
0032 pressure(1:6:m_t.gridset.gsize)=md.pressure;
0033
0034 if strcmpi(analysis_type,'thermalsteady'),
0035
0036
0037 disp(' computing temperature...');
0038 gridset=m_t.gridset;
0039 inputs=struct('pressure',pressure,'velocity',velocity,'dt',md.dt);
0040 [t_g loads_t melting_offset]=icethermal_core(m_t,analysis_type,inputs);
0041
0042
0043 disp(' computing melting...');
0044 gridset=m_m.gridset;
0045 inputs=struct('pressure',pressure,'temperature',t_g,'melting_offset',melting_offset);
0046 melting_g=icemelting_core(m_m,analysis_type,inputs);
0047
0048
0049 indx=1:6:gridset.gsize;
0050 indx=indx(m_t.tpart);
0051 md.temperature=t_g(indx);
0052 md.melting=melting_g(indx)*md.yts;
0053 else
0054
0055 nsteps=md.ndt/md.dt;
0056
0057
0058 temperature=zeros(gridset.gsize,1);
0059 melting=zeros(gridset.gsize,1);
0060 temperature(1:6:gridset.gsize)=md.temperature;
0061 melting(1:6:gridset.gsize)=md.melting/md.yts;
0062 soln.t_g=temperature;
0063 soln.melting_g=melting;
0064 soln.time=0;
0065
0066
0067 for n=1:nsteps,
0068
0069 soln(n+1).time=(n)*md.dt;
0070 disp(sprintf('\n%s%i/%i\n','time step: ',n,md.ndt/md.dt));
0071
0072
0073 disp(' computing temperature...');
0074 gridset=m_t.gridset;
0075 inputs=struct('pressure',pressure,'temperature',soln(n).t_g,'velocity',velocity,'dt',md.dt);
0076 [soln(n+1).t_g loads_t melting_offset]=icethermal_core(m_t,analysis_type,inputs);
0077
0078
0079 disp(' computing melting...');
0080 gridset=m_m.gridset;
0081 inputs=struct('pressure',pressure,'temperature',soln(n).t_g,'melting_offset',melting_offset,'dt',md.dt);
0082 soln(n+1).melting_g=icemelting_core(m_m,analysis_type,inputs);
0083
0084 end
0085
0086
0087 indx=1:6:m_t.gridset.gsize;
0088 indx=indx(m_t.tpart);
0089 solution=struct('time',{},'temperature',{},'melting',{});
0090
0091 for n=1:nsteps+1,
0092 solution(n).temperature=soln(n).t_g(indx);
0093 solution(n).melting=soln(n).melting_g(indx)*md.yts;
0094 solution(n).time=soln(n).time/md.yts;
0095 end
0096 md.thermaltransient_results=solution;
0097 end