0001 function [elements,grids,loads,constraints,materials,part,tpart]=ModelProcessorThermal(md,solutiontype);
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013 global cluster
0014
0015 if cluster,
0016
0017 element_partitioning=MeshPartition(md,numlabs);
0018 else
0019
0020 element_partitioning=ones(md.numberofelements,1);
0021 labindex=1;
0022 end
0023
0024
0025 elements=struct('element',cell(md.numberofelements,1));
0026 materials=struct('material',cell(md.numberofelements+1,1));
0027 mygrids=zeros(md.numberofgrids,1);
0028
0029
0030 pos=find(element_partitioning==labindex);
0031 [elements(pos).element]=deal(pentaelem);
0032
0033 elements(pos)=SetStructureField(elements(pos),'element','type','pentaelem');
0034 elements(pos)=SetStructureField(elements(pos),'element','id',pos);
0035 elements(pos)=SetStructureField(elements(pos),'element','g',md.elements(pos,1:6));
0036 elements(pos)=SetStructureField(elements(pos),'element','h',md.thickness(md.elements(pos,1:6)));
0037 elements(pos)=SetStructureField(elements(pos),'element','s',md.surface(md.elements(pos,1:6)));
0038 elements(pos)=SetStructureField(elements(pos),'element','b',md.bed(md.elements(pos,1:6)));
0039 elements(pos)=SetStructureField(elements(pos),'element','friction_type',md.drag_type);
0040 elements(pos)=SetStructureField(elements(pos),'element','k',md.drag(md.elements(pos,1:6)));
0041 elements(pos)=SetStructureField(elements(pos),'element','p',md.p(pos));
0042 elements(pos)=SetStructureField(elements(pos),'element','q',md.q(pos));
0043 elements(pos)=SetStructureField(elements(pos),'element','shelf',md.elementoniceshelf(pos));
0044 elements(pos)=SetStructureField(elements(pos),'element','onbed',md.elementonbed(pos));
0045 elements(pos)=SetStructureField(elements(pos),'element','onsurface',md.elementonsurface(pos));
0046 elements(pos)=SetStructureField(elements(pos),'element','meanvel',md.meanvel);
0047 elements(pos)=SetStructureField(elements(pos),'element','epsvel',md.epsvel);
0048 elements(pos)=SetStructureField(elements(pos),'element','acceleration',0);
0049 elements(pos)=SetStructureField(elements(pos),'element','collapse',0);
0050 elements(pos)=SetStructureField(elements(pos),'element','matid',pos);
0051 elements(pos)=SetStructureField(elements(pos),'element','geothermalflux',md.geothermalflux(md.elements(pos,1:6)));
0052 if strcmpi(solutiontype,'thermalsteady'),
0053 elements(pos)=SetStructureField(elements(pos),'element','thermal_steadystate',1);
0054 else
0055 elements(pos)=SetStructureField(elements(pos),'element','thermal_steadystate',0);
0056 end
0057
0058 [materials(pos).material]=deal(matice);
0059
0060 materials(pos)=SetStructureField(materials(pos),'material','id',pos);
0061 materials(pos)=SetStructureField(materials(pos),'material','B',md.B(md.elements(pos,1:6))*[1;1;1;1;1;1]/6);
0062 materials(pos)=SetStructureField(materials(pos),'material','n',md.n(pos));
0063
0064
0065 [materials(end).constants]=deal(matpar);
0066 materials(end)=SetStructureField(materials(end),'constants','g',md.g);
0067 materials(end)=SetStructureField(materials(end),'constants','rho_ice',md.rho_ice);
0068 materials(end)=SetStructureField(materials(end),'constants','rho_water',md.rho_water);
0069 materials(end)=SetStructureField(materials(end),'constants','thermalconductivity',md.thermalconductivity);
0070 materials(end)=SetStructureField(materials(end),'constants','heatcapacity',md.heatcapacity);
0071 materials(end)=SetStructureField(materials(end),'constants','latentheat',md.latentheat);
0072 materials(end)=SetStructureField(materials(end),'constants','beta',md.beta);
0073 materials(end)=SetStructureField(materials(end),'constants','meltingpoint',md.meltingpoint);
0074 materials(end)=SetStructureField(materials(end),'constants','mixed_layer_capacity',md.mixed_layer_capacity);
0075 materials(end)=SetStructureField(materials(end),'constants','thermal_exchange_velocity',md.thermal_exchange_velocity);
0076
0077 if cluster,
0078
0079 mygrids(md.elements(el3pos,:))=1;
0080 mygrids(md.elements(el6pos,:))=1;
0081 end
0082
0083 if cluster,
0084
0085
0086 bordergrids=double(gplus(mygrids)>1);
0087 else
0088 bordergrids=zeros(md.numberofgrids,1);
0089 end
0090
0091
0092 grids=struct('grid',cell(md.numberofgrids,1));
0093
0094 pos=[1:md.numberofgrids]';
0095 [grids(pos).grid]=deal(node);
0096 grids(pos)=SetStructureField(grids(pos),'grid','id',pos);
0097 grids(pos)=SetStructureField(grids(pos),'grid','x',md.x(pos));
0098 grids(pos)=SetStructureField(grids(pos),'grid','y',md.y(pos));
0099 grids(pos)=SetStructureField(grids(pos),'grid','z',md.z(pos));
0100 grids(pos)=SetStructureField(grids(pos),'grid','onbed',md.gridonbed(pos));
0101 grids(pos)=SetStructureField(grids(pos),'grid','border',bordergrids(pos));
0102
0103
0104
0105
0106
0107
0108
0109
0110
0111
0112
0113
0114 grids(pos)=SetStructureField(grids(pos),'grid','gridset','23456');
0115
0116
0117
0118
0119 loads=struct('load',cell(0,1));
0120
0121
0122 count=1;
0123 for i=1:md.numberofgrids,
0124 if ~md.gridondirichlet_thermal(i),
0125 pengridobject=pengrid;
0126 pengridobject.id=count;
0127 pengridobject.dof=1;
0128 pengridobject.grid=i;
0129 pengridobject.active=0;
0130 if strcmpi(solutiontype,'thermalsteady'),
0131 pengridobject.thermal_steadystate=1;
0132 else
0133 pengridobject.thermal_steadystate=0;
0134 end
0135 pengridobject.penalty_offset=md.penalty_offset;
0136
0137 loads(count).load=pengridobject;
0138 count=count+1;
0139 end
0140 end
0141
0142
0143 spcs=find(md.gridondirichlet_thermal);
0144 constraints=struct('constraint',cell(length(spcs),1));
0145
0146 count=1;
0147 for i=1:md.numberofgrids,
0148 if md.gridondirichlet_thermal(i),
0149 constraints(count).constraint=spc;
0150 constraints(count).constraint.grid=i;
0151 constraints(count).constraint.dof=1;
0152 constraints(count).constraint.value=md.dirichletvalues_thermal(i,1);
0153 count=count+1;
0154 end
0155 end
0156
0157
0158
0159 [part,tpart]=PartitioningVector(md,grids);
0160
0161
0162 end