Index: /issm/trunk/test/Validation/ThermalTests/Melting/Square.par =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Melting/Square.par (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Melting/Square.par (revision 4876) @@ -3,4 +3,6 @@ %parallelization md.cluster='none'; +md.np=3; +md.dt=0; disp(' creating thickness'); @@ -39,9 +41,10 @@ disp(' boundary conditions for diagnostic model'); -md=SetIceShelfBC(md,'Front.exp'); +md=SetMarineIceSheetBC(md,'Front.exp'); disp(' boundary conditions for thermal model'); -md.gridondirichlet_thermal=ones(md.numberofgrids,1); %surface temperature -md.dirichletvalues_thermal=md.observed_temperature; +md.temperature=md.observed_temperature; +md.spctemperature(:,1)=ones(md.numberofgrids,1); %surface temperature +md.spctemperature(:,2)=md.observed_temperature; md.geothermalflux=zeros(md.numberofgrids,1); pos=find(md.elementonicesheet);md.geothermalflux(md.elements(pos,:))=1*10^-3; %1 mW/m^2 Index: /issm/trunk/test/Validation/ThermalTests/Melting/runme.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Melting/runme.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Melting/runme.m (revision 4876) @@ -22,9 +22,10 @@ %modeled results -%md=solve(md,'analysis_type','thermal','sub_analysis_type','steady'); -md=solve(md,'analysis_type','thermal','sub_analysis_type','steady'); +md=solve(md,'analysis_type',ThermalSolutionEnum); %plot results -relative=abs((md.results.thermal.melting-melting)./melting)*100; -relative(find(md.results.thermal.melting==melting))=0; -plotmodel(md,'data','melting','title','modeled melting','data',melting,'title','analytical melting','data',md.results.thermal.melting-melting,'title','absolute error','data',relative,'title','relative error [%]','layer#all',1) +comp_melting=zeros(md.numberofgrids,1); +comp_melting(md.results.ThermalSolution.MeltingRate.index)=md.results.ThermalSolution.MeltingRate.value; +relative=abs((comp_melting-melting)./melting)*100; +relative(find(comp_melting==melting))=0; +plotmodel(md,'data',comp_melting,'title','modeled melting','data',melting,'title','analytical melting','data',comp_melting-melting,'title','absolute error','data',relative,'title','relative error [%]','layer#all',1) Index: /issm/trunk/test/Validation/ThermalTests/Melting/thermalboundarycondition.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Melting/thermalboundarycondition.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Melting/thermalboundarycondition.m (revision 4876) @@ -2,5 +2,5 @@ %Ok, start defining model BC disp(' boundary conditions for thermal model '); -if strcmp(md.type,'3d') +if md.dim==3, % pos2=find(md.elementonsurface); md.gridondirichlet_thermal(md.elements(pos2,4:6))=1; md.dirichletvalues_thermal(md.elements(pos2,4:6))=0; md.vx=zeros(md.numberofgrids,1); Index: /issm/trunk/test/Validation/ThermalTests/Simpleadvection/Square.par =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simpleadvection/Square.par (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simpleadvection/Square.par (revision 4876) @@ -3,4 +3,5 @@ %parallelization md.cluster='none'; +md.dt=0; disp(' creating thickness'); @@ -37,11 +38,14 @@ %Deal with boundary conditions: - disp(' boundary conditions for diagnostic model'); -md=SetIceShelfBC(md,'Front.exp'); +md=SetMarineIceSheetBC(md,'Front.exp'); disp(' boundary conditions for thermal model'); -md.gridondirichlet_thermal=ones(md.numberofgrids,1); %surface temperature -md.dirichletvalues_thermal=md.observed_temperature; +md.temperature=md.observed_temperature; +md.spctemperature(:,1)=ones(md.numberofgrids,1); %surface temperature +md.spctemperature(:,2)=md.observed_temperature; md.geothermalflux=zeros(md.numberofgrids,1); pos=find(md.elementonicesheet);md.geothermalflux(md.elements(pos,:))=50*10^-3; %50 mW/m^2 + +disp(' cluster properties'); +md.np=3; Index: /issm/trunk/test/Validation/ThermalTests/Simpleadvection/runme.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simpleadvection/runme.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simpleadvection/runme.m (revision 4876) @@ -23,8 +23,10 @@ %modeled results -md=solve(md,'analysis_type','thermal','sub_analysis_type','steady','package','ice'); +md=solve(md,'analysis_type',ThermalSolutionEnum); %plot results -relative=abs((md.results.thermal.temperature-md.observed_temperature)./md.observed_temperature)*100; -relative(find(md.results.thermal.temperature==md.observed_temperature))=0; -plotmodel(md,'data',md.results.thermal.temperature,'title','modeled temperature','data','observed_temperature','view',3,'title','analytical temperature','view',3,'data',md.results.thermal.temperature-md.observed_temperature,'title','absolute error','view',3,'data',relative,'title','relative error [%]','view',3) +comp_temp=zeros(md.numberofgrids,1); +comp_temp(md.results.ThermalSolution.Temperature.index)=md.results.ThermalSolution.Temperature.value; +relative=abs((comp_temp-md.observed_temperature)./md.observed_temperature)*100; +relative(find(comp_temp==md.observed_temperature))=0; +plotmodel(md,'data',comp_temp,'title','modeled temperature','data','observed_temperature','view',3,'title','analytical temperature','view',3,'data',comp_temp-md.observed_temperature,'title','absolute error','view',3,'data',relative,'title','relative error [%]','view',3) Index: /issm/trunk/test/Validation/ThermalTests/Simpleadvection/thermalboundarycondition.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simpleadvection/thermalboundarycondition.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simpleadvection/thermalboundarycondition.m (revision 4876) @@ -2,7 +2,7 @@ %Ok, start defining model BC disp(' boundary conditions for thermal model '); -if strcmp(md.type,'3d') - pos1=find(md.elementonbed); md.gridondirichlet_thermal(md.elements(pos1,1:3))=1; md.dirichletvalues_thermal(md.elements(pos1,1:3))=10; - pos2=find(md.elementonsurface); md.gridondirichlet_thermal(md.elements(pos2,4:6))=1; md.dirichletvalues_thermal(md.elements(pos2,4:6))=0; +if md.dim==3, + pos1=find(md.elementonbed); md.spctemperature(md.elements(pos1,1:3),1)=1; md.spctemperature(md.elements(pos1,1:3),2)=10; + pos2=find(md.elementonsurface); md.spctemperature(md.elements(pos2,4:6),1)=1; md.spctemperature(md.elements(pos2,4:6),2)=0; md.vx=zeros(md.numberofgrids,1); md.vy=zeros(md.numberofgrids,1); Index: /issm/trunk/test/Validation/ThermalTests/Simpleconduction/Square.par =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simpleconduction/Square.par (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simpleconduction/Square.par (revision 4876) @@ -1,47 +1,49 @@ %Ok, start defining model parameters here - %parallelization - md.cluster='none'; +md.dt=0; - disp(' creating thickness'); - h=1000; - md.thickness=h*ones(md.numberofgrids,1); - md.firn_layer=10*ones(md.numberofgrids,1); - md.bed=zeros(md.numberofgrids,1); - md.surface=md.thickness; +%parallelization +md.np=3; - disp(' creating velocities'); - md.vx_obs=zeros(md.numberofgrids,1); - md.vy_obs=zeros(md.numberofgrids,1); - md.vel_obs=sqrt(md.vx_obs.^2+md.vy_obs.^2); - - disp(' creating drag'); - md.drag_type=2; %0 none 1 plastic 2 viscous - md.drag_coefficient=200*ones(md.numberofgrids,1); %q=1. - %Take care of iceshelves: no basal drag - pos=find(md.elementoniceshelf); - md.drag_coefficient(md.elements(pos,:))=0; - md.drag_p=ones(md.numberofelements,1); - md.drag_q=ones(md.numberofelements,1); +disp(' creating thickness'); +h=1000; +md.thickness=h*ones(md.numberofgrids,1); +md.firn_layer=10*ones(md.numberofgrids,1); +md.bed=zeros(md.numberofgrids,1); +md.surface=md.thickness; - disp(' creating temperatures'); - md.observed_temperature=(273-20)*ones(md.numberofgrids,1); +disp(' creating velocities'); +md.vx_obs=zeros(md.numberofgrids,1); +md.vy_obs=zeros(md.numberofgrids,1); +md.vel_obs=sqrt(md.vx_obs.^2+md.vy_obs.^2); - disp(' creating flow law paramter'); - md.rheology_B=paterson(md.observed_temperature); - md.rheology_n=3*ones(md.numberofelements,1); - - disp(' creating accumulation rates'); - md.accumulation_rate=ones(md.numberofgrids,1)/md.yts; %1m/a - md.melting_rate=0*ones(md.numberofgrids,1)/md.yts; %1m/a +disp(' creating drag'); +md.drag_type=2; %0 none 1 plastic 2 viscous +md.drag_coefficient=200*ones(md.numberofgrids,1); %q=1. +%Take care of iceshelves: no basal drag +pos=find(md.elementoniceshelf); +md.drag_coefficient(md.elements(pos,:))=0; +md.drag_p=ones(md.numberofelements,1); +md.drag_q=ones(md.numberofelements,1); - %Deal with boundary conditions: - - disp(' boundary conditions for diagnostic model'); - md=SetIceShelfBC(md,'Front.exp'); - - disp(' boundary conditions for thermal model'); - md.gridondirichlet_thermal=ones(md.numberofgrids,1); %surface temperature - md.dirichletvalues_thermal=md.observed_temperature; - md.geothermalflux=zeros(md.numberofgrids,1); - pos=find(md.elementonicesheet);md.geothermalflux(md.elements(pos,:))=50*10^-3; %50 mW/m^2 +disp(' creating temperatures'); +md.observed_temperature=(273-20)*ones(md.numberofgrids,1); + +disp(' creating flow law paramter'); +md.rheology_B=paterson(md.observed_temperature); +md.rheology_n=3*ones(md.numberofelements,1); + +disp(' creating accumulation rates'); +md.accumulation_rate=ones(md.numberofgrids,1)/md.yts; %1m/a +md.melting_rate=0*ones(md.numberofgrids,1)/md.yts; %1m/a + +%Deal with boundary conditions: + +disp(' boundary conditions for diagnostic model'); +md=SetMarineIceSheetBC(md,'Front.exp'); + +disp(' boundary conditions for thermal model'); +md.spctemperature(:,1)=ones(md.numberofgrids,1); %surface temperature +md.spctemperature(:,2)=md.observed_temperature; +md.geothermalflux=zeros(md.numberofgrids,1); +pos=find(md.elementonicesheet);md.geothermalflux(md.elements(pos,:))=50*10^-3; %50 mW/m^2 Index: /issm/trunk/test/Validation/ThermalTests/Simpleconduction/runme.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simpleconduction/runme.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simpleconduction/runme.m (revision 4876) @@ -22,8 +22,10 @@ %modeled results -md=solve(md,'analysis_type','thermal','sub_analysis_type','steady'); +md=solve(md,'analysis_type',ThermalSolutionEnum); %plot results -relative=abs((md.results.thermal.temperature-md.observed_temperature)./md.observed_temperature)*100; -relative(find(md.results.thermal.temperature==md.observed_temperature))=0; -plotmodel(md,'data','temperature','title','modeled temperature','data','observed_temperature','view',3,'title','analytical temperature','view',3,'data',md.results.thermal.temperature-md.observed_temperature,'title','absolute error','view',3,'data',relative,'title','relative error [%]','view',3) +comp_temp=zeros(md.numberofgrids,1); +comp_temp(md.results.ThermalSolution.Temperature.index)=md.results.ThermalSolution.Temperature.value; +relative=abs((comp_temp-md.observed_temperature)./md.observed_temperature)*100; +relative(find(comp_temp==md.observed_temperature))=0; +plotmodel(md,'data',comp_temp,'title','modeled temperature','data','observed_temperature','view',3,'title','analytical temperature','view',3,'data',comp_temp-md.observed_temperature,'title','absolute error','view',3,'data',relative,'title','relative error [%]','view',3) Index: /issm/trunk/test/Validation/ThermalTests/Simpleconduction/thermalboundarycondition.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simpleconduction/thermalboundarycondition.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simpleconduction/thermalboundarycondition.m (revision 4876) @@ -2,7 +2,7 @@ %Ok, start defining model BC disp(' boundary conditions for thermal model '); -if strcmp(md.type,'3d') - pos1=find(md.elementonbed); md.gridondirichlet_thermal(md.elements(pos1,1:3))=1; md.dirichletvalues_thermal(md.elements(pos1,1:3))=10; - pos2=find(md.elementonsurface); md.gridondirichlet_thermal(md.elements(pos2,4:6))=1; md.dirichletvalues_thermal(md.elements(pos2,4:6))=0; +if md.dim==3, + pos1=find(md.elementonbed); md.spctemperature(md.elements(pos1,1:3),1)=1; md.spctemperature(md.elements(pos1,1:3),2)=10; + pos2=find(md.elementonsurface); md.spctemperature(md.elements(pos2,4:6),1)=1; md.spctemperature(md.elements(pos2,4:6),2)=0; md.vx=zeros(md.numberofgrids,1); md.vy=zeros(md.numberofgrids,1); Index: /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/Square.par =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/Square.par (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/Square.par (revision 4876) @@ -1,6 +1,7 @@ %Ok, start defining model parameters here +md.dt=0; %parallelization -md.cluster='none'; +md.np=3; disp(' creating thickness'); @@ -39,9 +40,9 @@ disp(' boundary conditions for diagnostic model'); -md=SetIceShelfBC(md,'Front.exp'); +md=SetMarineIceSheetBC(md,'Front.exp'); disp(' boundary conditions for thermal model'); -md.gridondirichlet_thermal=ones(md.numberofgrids,1); %surface temperature -md.dirichletvalues_thermal=md.observed_temperature; +md.spctemperature(:,1)=ones(md.numberofgrids,1); %surface temperature +md.spctemperature(:,2)=md.observed_temperature; md.geothermalflux=zeros(md.numberofgrids,1); pos=find(md.elementonicesheet);md.geothermalflux(md.elements(pos,:))=0.1; %100 mW/m^2 Index: /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/runme.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/runme.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/runme.m (revision 4876) @@ -22,8 +22,10 @@ %modeled results -md=solve(md,'analysis_type','thermal','sub_analysis_type','steady'); +md=solve(md,'analysis_type',ThermalSolutionEnum); %plot results -relative=abs((md.results.thermal.temperature-md.observed_temperature)./md.observed_temperature)*100; -relative(find(md.results.thermal.temperature==md.observed_temperature))=0; -plotmodel(md,'data','temperature','title','modeled temperature','data','observed_temperature','view',3,'title','analytical temperature','view',3,'data',md.results.thermal.temperature-md.observed_temperature,'title','absolute error','view',3,'data',relative,'title','relative error [%]','view',3) +comp_temp=zeros(md.numberofgrids,1); +comp_temp(md.results.ThermalSolution.Temperature.index)=md.results.ThermalSolution.Temperature.value; +relative=abs((comp_temp-md.observed_temperature)./md.observed_temperature)*100; +relative(find(comp_temp==md.observed_temperature))=0; +plotmodel(md,'data',comp_temp,'title','modeled temperature','data','observed_temperature','view',3,'title','analytical temperature','view',3,'data',comp_temp-md.observed_temperature,'title','absolute error','view',3,'data',relative,'title','relative error [%]','view',3) Index: /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/thermalboundarycondition.m =================================================================== --- /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/thermalboundarycondition.m (revision 4875) +++ /issm/trunk/test/Validation/ThermalTests/Simplegeothermalflux/thermalboundarycondition.m (revision 4876) @@ -2,6 +2,6 @@ %Ok, start defining model BC disp(' boundary conditions for thermal model '); -if strcmp(md.type,'3d') - pos2=find(md.elementonsurface); md.gridondirichlet_thermal(md.elements(pos2,4:6))=1; md.dirichletvalues_thermal(md.elements(pos2,4:6))=0; +if md.dim==3, + pos2=find(md.elementonsurface); md.spctemperature(md.elements(pos2,4:6),1)=1; md.spctemperature(md.elements(pos2,4:6),2)=0; md.vx=zeros(md.numberofgrids,1); md.vy=zeros(md.numberofgrids,1);