cielothermal

PURPOSE ^

thermal solution: steady state and transient

SYNOPSIS ^

function md=cielothermal(md,solutiontype)

DESCRIPTION ^

thermal  solution: steady state and transient

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function md=cielothermal(md,solutiontype)
0002 %thermal  solution: steady state and transient
0003 
0004     analysis_t=solutiontype;
0005     analysis_m='melting';
0006     
0007     %configure, and initialise data model and parameter defaults:
0008     SetParameterDefaults;
0009 
0010     %Build thermal and melting models
0011     m_t=CreateFEMModel(md,analysis_t);
0012     m_m=CreateFEMModel(md,analysis_m);
0013 
0014 
0015     %Build partitioning vectors.
0016     indx=1:6:mt.uset.gsize; indx=indx(tpart);
0017 
0018     % figure out number of dof: just for information purposes.
0019     md.dof=m_t.uset.fsize; %biggest dof number
0020 
0021     %initialize velocity,pressure and material properties
0022     velocity=zeros(m_t.uset.gsize,1);
0023     velocity(1:6:m_t.uset.gsize)=md.vx(part)/md.yts;
0024     velocity(2:6:m_t.uset.gsize)=md.vy(part)/md.yts;
0025     velocity(3:6:m_t.uset.gsize)=md.vz(part)/md.yts;
0026 
0027     pressure=zeros(m_t.uset.gsize,1);
0028     pressure(1:6:m_t.uset.gsize)=md.pressure(part);
0029 
0030     thermalparams.beta=md.beta;
0031     thermalparams.meltingpoint=md.meltingpoint;
0032     thermalparams.latentheat=md.latentheat;
0033     thermalparams.heatcapacity=md.heatcapacity;
0034     thermalparams.penalty=md.penalty_thermal;
0035     thermalparams.min_thermal_constraints=md.min_thermal_constraints;
0036     thermalparams.mintemp=min(md.dirichletvalues_thermal(find(md.gridondirichlet_thermal)));
0037 
0038     if strcmpi(analysis_t,'thermalsteady'),
0039 
0040         inputs=struct('pressure',pressure,'velocity',velocity,'dt',md.dt);
0041 
0042         %Call core thermal computation
0043         disp(sprintf('\n%s\n','computing temperature...'));
0044         [t_g,m_t]=cielothermal_core(m_t,m_t.params,thermalparams,inputs,analysis_t);
0045 
0046         %Call core melting computation:
0047         %disp(sprintf('\n%s\n','computing melting...'));
0048         %SetUset(m_m);
0049         %inputs=struct('pressure',pressure,'temperature',rt_g,'dt',md.dt);
0050         %rmelting_g=cielomelting_core(m_m,m_m.params,thermalparams,inputs,analysis_m);
0051 
0052         %load results onto model
0053         md.temperature=t_g(indx);
0054         %melting_g=IMdb('select matrix from rmelting_g');
0055         %md.melting=melting_g(indx)*md.yts; %from m/s to m/a
0056     else
0057 
0058         nsteps=md.ndt/md.dt;
0059 
0060         %initialize temperature and melting
0061         temperature=zeros(m_t.uset.gsize,1);
0062         temperature(1:6:m_t.uset.gsize)=md.temperature(part);
0063         soln.t_g=temperature;
0064         %melting=zeros(uset.gsize,1);
0065         %melting(1:6:uset.gsize)=md.melting(part)/md.yts; %in m/s
0066         %soln.melting_g=melting;
0067 
0068         for n=1:nsteps, 
0069 
0070             disp(sprintf('\n%s%i/%i\n','time step: ',n,md.ndt/md.dt));
0071             %Call core thermal computation
0072             disp('   computing temperature...');
0073             inputs=struct('pressure',pressure,'temperature',soln(n).t_g,'velocity',velocity,'dt',md.dt);
0074             [soln(n+1).t_g,m_t]=cielothermal_core(m_t,m_t.params,thermalparams,inputs,analysis_t);
0075             
0076             %Call core melting computation
0077             %disp('   computing melting...');
0078             %SetUset(m_m);
0079             %inputs=struct('pressure',pressure,'temperature',soln(n).t_g,'dt',md.dt);
0080             %rmelting_g=cielomelting_core(m_m,m_m.params,thermalparams,inputs,analysis_m);
0081             %soln(n+1).melting_g=IMdb('select matrix from rmelting_g');
0082             
0083         end
0084         
0085         %Wrap up
0086         solution_temperature=struct('temperature',{});
0087         for n=1:nsteps+1,
0088             solution_temperature(n).temperature=soln(n).t_g(indx);
0089         end
0090         md.temperature=solution_temperature;
0091         
0092         %solution_melting=struct('melting',{});
0093         %for n=1:nsteps+1,
0094         %    solution_melting(n).melting=soln(n).melting_g(indx)*md.yts; %in m/a
0095         %end
0096         %md.melting=solution_melting;
0097     end
0098

Generated on Sun 29-Mar-2009 20:22:55 by m2html © 2003