Index: /issm/trunk/src/m/solutions/ice/ModelProcessorDiagnosticHoriz.m =================================================================== --- /issm/trunk/src/m/solutions/ice/ModelProcessorDiagnosticHoriz.m (revision 1751) +++ /issm/trunk/src/m/solutions/ice/ModelProcessorDiagnosticHoriz.m (revision 1752) @@ -159,43 +159,43 @@ %icefront -loads=struct('load',cell(length(md.segmentonneumann_diag),1)); - -for i=1:size(md.segmentonneumann_diag,1), +loads=struct('load',cell(length(md.pressureload),1)); + +for i=1:size(md.pressureload,1), if strcmpi(md.type,'3d'), - if md.elements_type(md.segmentonneumann_diag(i,end))==MacAyealFormulationEnum(), + if md.elements_type(md.pressureload(i,end))==MacAyealFormulationEnum(), loads(i).load=icefront; - loads(i).load.eid=find(el6pos==md.segmentonneumann_diag(i,end)); - loads(i).load.g(1)=md.segmentonneumann_diag(i,1); - loads(i).load.g(2)=md.segmentonneumann_diag(i,2); + loads(i).load.eid=find(el6pos==md.pressureload(i,end)); + loads(i).load.g(1)=md.pressureload(i,1); + loads(i).load.g(2)=md.pressureload(i,2); loads(i).load.rho_water=md.rho_water; loads(i).load.type='segment'; - elseif md.elements_type(md.segmentonneumann_diag(i,end))==PattynFormulationEnum(), + elseif md.elements_type(md.pressureload(i,end))==PattynFormulationEnum(), %build a quad ice front for the penta element loads(i).load=icefront; - loads(i).load.eid=find(el6pos==md.segmentonneumann_diag(i,end)); - loads(i).load.g(1)=md.segmentonneumann_diag(i,1); - loads(i).load.g(2)=md.segmentonneumann_diag(i,2); - loads(i).load.g(3)=md.segmentonneumann_diag(i,3); - loads(i).load.g(4)=md.segmentonneumann_diag(i,4); + loads(i).load.eid=find(el6pos==md.pressureload(i,end)); + loads(i).load.g(1)=md.pressureload(i,1); + loads(i).load.g(2)=md.pressureload(i,2); + loads(i).load.g(3)=md.pressureload(i,3); + loads(i).load.g(4)=md.pressureload(i,4); loads(i).load.type='quad'; else - if ~(md.elements_type(md.segmentonneumann_diag(i,end))==HutterFormulationEnum()), + if ~(md.elements_type(md.pressureload(i,end))==HutterFormulationEnum()), error('diagnostic error message: unsupported element type'); end end else - if md.elements_type(md.segmentonneumann_diag(i,end))==MacAyealFormulationEnum(), + if md.elements_type(md.pressureload(i,end))==MacAyealFormulationEnum(), loads(i).load=icefront; - loads(i).load.eid=find(el3pos==md.segmentonneumann_diag(i,end)); - loads(i).load.g(1)=md.segmentonneumann_diag(i,1); - loads(i).load.g(2)=md.segmentonneumann_diag(i,2); + loads(i).load.eid=find(el3pos==md.pressureload(i,end)); + loads(i).load.g(1)=md.pressureload(i,1); + loads(i).load.g(2)=md.pressureload(i,2); loads(i).load.rho_water=md.rho_water; loads(i).load.type='segment'; else - if ~(md.elements_type(md.segmentonneumann_diag(i,end))==HutterFormulationEnum()), + if ~(md.elements_type(md.pressureload(i,end))==HutterFormulationEnum()), error('diagnostic error message: unsupported element type'); end @@ -207,5 +207,5 @@ %Initialize constraints structure -constraints=struct('constraint',cell(2*length(find(md.gridondirichlet_diag)),1)); +constraints=struct('constraint',cell(length(find(md.spcvelocity(:,1)))+length(find(md.spcvelocity(:,2))),1)); count=1; @@ -228,19 +228,22 @@ count=count+1; - elseif (md.gridonmacayeal(i) | md.gridonpattyn(i)) & md.gridondirichlet_diag(i), - - %constrain first dof - constraints(count).constraint=spc; - constraints(count).constraint.grid=i; - constraints(count).constraint.dof=1; - constraints(count).constraint.value=md.dirichletvalues_diag(i,1)/md.yts; %in m/s - count=count+1; - - %constrain second dof - constraints(count).constraint=spc; - constraints(count).constraint.grid=i; - constraints(count).constraint.dof=2; - constraints(count).constraint.value=md.dirichletvalues_diag(i,2)/md.yts; %in m/s - count=count+1; + else + if (md.gridonmacayeal(i) | md.gridonpattyn(i)) & md.spcvelocity(i,1), + %constrain first dof + constraints(count).constraint=spc; + constraints(count).constraint.grid=i; + constraints(count).constraint.dof=1; + constraints(count).constraint.value=md.spcvelocity(i,4)/md.yts; %in m/s + count=count+1; + end + + if (md.gridonmacayeal(i) | md.gridonpattyn(i)) & md.spcvelocity(i,2), + %constrain second dof + constraints(count).constraint=spc; + constraints(count).constraint.grid=i; + constraints(count).constraint.dof=2; + constraints(count).constraint.value=md.spcvelocity(i,5)/md.yts; %in m/s + count=count+1; + end end end Index: /issm/trunk/src/m/solutions/ice/ModelProcessorDiagnosticStokes.m =================================================================== --- /issm/trunk/src/m/solutions/ice/ModelProcessorDiagnosticStokes.m (revision 1751) +++ /issm/trunk/src/m/solutions/ice/ModelProcessorDiagnosticStokes.m (revision 1752) @@ -107,21 +107,21 @@ %icefront -if isnan(md.segmentonneumann_diag_stokes), - segmentonneumann_diag_stokes=[]; +if isnan(md.pressureload_stokes), + pressureload_stokes=[]; else - segmentonneumann_diag_stokes=md.segmentonneumann_diag_stokes; + pressureload_stokes=md.pressureload_stokes; end -length_segmentonneumann_diag_stokes= size(segmentonneumann_diag_stokes,1); -loads=struct('load',cell(length_segmentonneumann_diag_stokes,1)); +length_pressureload_stokes= size(pressureload_stokes,1); +loads=struct('load',cell(length_pressureload_stokes,1)); -for i=1:length_segmentonneumann_diag_stokes, +for i=1:length_pressureload_stokes, %build a quad ice front for the penta element loads(i).load=icefront; - loads(i).load.eid=find(segmentonneumann_diag_stokes(i,end)==find(md.elements_type(:,2)==StokesFormulationEnum())); %elements contain only stokes elements, so we have to renumbered - loads(i).load.g(1)=segmentonneumann_diag_stokes(i,1); - loads(i).load.g(2)=segmentonneumann_diag_stokes(i,2); - loads(i).load.g(3)=segmentonneumann_diag_stokes(i,3); - loads(i).load.g(4)=segmentonneumann_diag_stokes(i,4); + loads(i).load.eid=find(pressureload_stokes(i,end)==find(md.elements_type(:,2)==StokesFormulationEnum())); %elements contain only stokes elements, so we have to renumbered + loads(i).load.g(1)=pressureload_stokes(i,1); + loads(i).load.g(2)=pressureload_stokes(i,2); + loads(i).load.g(3)=pressureload_stokes(i,3); + loads(i).load.g(4)=pressureload_stokes(i,4); loads(i).load.rho_water=md.rho_water; loads(i).load.type='quad'; @@ -144,5 +144,5 @@ %Single point constraints: -spcs=find(md.gridondirichlet_diag); +spcs=find(md.spcvelocity(:,1)); constraints=struct('constraint',cell(3*length(spcs),1)); Index: /issm/trunk/src/m/solutions/ice/ModelProcessorPrognostic.m =================================================================== --- /issm/trunk/src/m/solutions/ice/ModelProcessorPrognostic.m (revision 1751) +++ /issm/trunk/src/m/solutions/ice/ModelProcessorPrognostic.m (revision 1752) @@ -95,9 +95,9 @@ %Initialize constraints structure -constraints=struct('constraint',cell(length(find(md.gridondirichlet_prog)),1)); +constraints=struct('constraint',cell(length(find(md.spcthickness(:,1))),1)); -pos=find(md.gridondirichlet_prog); +pos=find(md.spcthickness(:,1)); count=1; -for i=1:length(find(md.gridondirichlet_prog)), +for i=1:length(find(md.spcthickness(:,1))), %constrain the thickness @@ -105,5 +105,5 @@ constraints(count).constraint.grid=pos(i); constraints(count).constraint.dof=1; - constraints(count).constraint.value=md.dirichletvalues_prog(pos(i)); %in m + constraints(count).constraint.value=md.spcthickness(pos(i),2); %in m count=count+1; end Index: /issm/trunk/src/m/solutions/ice/ModelProcessorThermal.m =================================================================== --- /issm/trunk/src/m/solutions/ice/ModelProcessorThermal.m (revision 1751) +++ /issm/trunk/src/m/solutions/ice/ModelProcessorThermal.m (revision 1752) @@ -102,5 +102,5 @@ count=1; for i=1:md.numberofgrids, - if ~md.gridondirichlet_thermal(i), %No penalty applied on spc grids! + if ~md.spctemperature(i,1), %No penalty applied on spc grids! pengridobject=pengrid; pengridobject.id=count; @@ -121,14 +121,14 @@ %Single point constraints: -spcs=find(md.gridondirichlet_thermal); +spcs=find(md.spctemperature(:,1)); constraints=struct('constraint',cell(length(spcs),1)); count=1; for i=1:md.numberofgrids, - if md.gridondirichlet_thermal(i), + if md.spctemperature(i,1), constraints(count).constraint=spc; constraints(count).constraint.grid=i; constraints(count).constraint.dof=1; - constraints(count).constraint.value=md.dirichletvalues_thermal(i,1); + constraints(count).constraint.value=md.spctemperature(i,2); count=count+1; end Index: /issm/trunk/src/m/solutions/ice/VelocityExtrude.m =================================================================== --- /issm/trunk/src/m/solutions/ice/VelocityExtrude.m (revision 1751) +++ /issm/trunk/src/m/solutions/ice/VelocityExtrude.m (revision 1752) @@ -13,5 +13,5 @@ %Find list of 2d grids that belong to the collapsed macayeal elements. %And remove grids on the border. - %These grids are on macayeal but not on pattyn nor hutter nor dirichlet + %These grids are on macayeal but not on pattyn nor hutter nor spc grids2d=md.gridonmacayeal(1:md.numberofgrids2d);% & ~(md.gridonhutter(1:md.numberofgrids2d) | md.gridonpattyn(1:md.numberofgrids2d) | md.gridondirichlet_diag(1:md.numberofgrids2d)); Index: /issm/trunk/src/m/solutions/ice/diagnostic3d.m =================================================================== --- /issm/trunk/src/m/solutions/ice/diagnostic3d.m (revision 1751) +++ /issm/trunk/src/m/solutions/ice/diagnostic3d.m (revision 1752) @@ -68,5 +68,5 @@ inputs.velocity_average=velocity_average; -%update dirichlet boundary conditions with these new base vertical velocities +%update spc boundary conditions with these new base vertical velocities m_dv=fem.m_dv; gridset=m_dv.gridset; @@ -105,5 +105,5 @@ inputs.velocity=u_g; - %update dirichlet boundary conditions with the velocities computed previously + %update spc boundary conditions with the velocities computed previously gridset=m_ds.gridset; m_ds.y_g=u_g; m_ds.ys=Reducevector_g(m_ds.y_g);