Index: /issm/trunk-jpl/src/m/classes/damage.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/damage.m	(revision 16160)
+++ /issm/trunk-jpl/src/m/classes/damage.m	(revision 16160)
@@ -0,0 +1,97 @@
+%DAMAGEICE class definition
+%
+%   Usage:
+%      damage=damage();
+
+classdef damage
+	properties (SetAccess=public)  
+	% {{{
+		%damage 
+		isdamage                   = NaN;
+		D                   = NaN;
+		law                 = '';
+		
+		%parameters for law 'initial': 
+		stress_threshold    = NaN;
+		c1                  = NaN;
+		c2                  = NaN;
+		c3                  = NaN;
+		c4                  = NaN;
+
+	end % }}}
+	methods
+		function obj = damage(varargin) % {{{
+			switch nargin
+				case 0
+					obj=setdefaultparameters(obj);
+				case 1
+					inputstruct=varargin{1};
+					list1 = properties('damage');
+					list2 = fieldnames(inputstruct);
+					for i=1:length(list1)
+						fieldname = list1{i};
+						if ismember(fieldname,list2),
+							obj.(fieldname) = inputstruct.(fieldname);
+						end
+					end
+				otherwise
+					error('constructor not supported');
+			end
+		end % }}}
+		function obj = setdefaultparameters(obj) % {{{
+
+			%damage parameters: 
+			obj.isdamage=0;
+			obj.D=0;
+			obj.law='undamaged';
+			
+			%damage parameters for 'initial' law of damage propagation
+			obj.stress_threshold=0;
+			obj.c1=0;
+			obj.c2=0;
+			obj.c3=0;
+			obj.c4=0;
+
+		end % }}}
+		function md = checkconsistency(obj,md,solution,analyses) % {{{
+			
+			md = checkfield(md,'materials.D','>=0',0,'size',[md.mesh.numberofvertices 1]);
+			md = checkfield(md,'materials.law','values',{'undamaged','pralong'});
+			if strcmpi(obj.law,'pralong'),
+				md = checkfield(md,'materials.c1','>=',0);
+				md = checkfield(md,'materials.c2','>=',0);
+				md = checkfield(md,'materials.c3','>=',0);
+				md = checkfield(md,'materials.c4','>=',0);
+				md = checkfield(md,'materials.stress_threshold','>=',0);
+			end
+
+		end % }}}
+		function disp(obj) % {{{
+			disp(sprintf('   Damage:\n'));
+
+			fielddisplay(obj,'D','damage tensor (scalar)');
+			fielddisplay(obj,'law','damage law (string) from {''undamaged'',''pralong''}');
+			if strcmpi(obj.law,'pralong'),
+				fielddisplay(obj,'c1','damage parameter 1');
+				fielddisplay(obj,'c2','damage parameter 2');
+				fielddisplay(obj,'c3','damage parameter 3');
+				fielddisplay(obj,'c4','damage parameter 4');
+				fielddisplay(obj,'stress_threshold','damage stress threshold [Pa]');
+			end
+
+		end % }}}
+		function marshall(obj,md,fid) % {{{
+			
+			WriteData(fid,'object',obj,'class','materials','fieldname','D','format','DoubleMat','mattype',1);
+			WriteData(fid,'object',obj,'class','materials','fieldname','law','format','String');
+			if strcmpi(obj.law,'pralong'),
+				WriteData(fid,'object',obj,'class','materials','fieldname','c1','format','Double');
+				WriteData(fid,'object',obj,'class','materials','fieldname','c2','format','Double');
+				WriteData(fid,'object',obj,'class','materials','fieldname','c3','format','Double');
+				WriteData(fid,'object',obj,'class','materials','fieldname','c4','format','Double');
+				WriteData(fid,'object',obj,'class','materials','fieldname','stress_threshold','format','Double');
+			end
+
+		end % }}}
+	end
+end
Index: /issm/trunk-jpl/src/m/classes/matdamageice.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/matdamageice.m	(revision 16159)
+++ /issm/trunk-jpl/src/m/classes/matdamageice.m	(revision 16160)
@@ -163,5 +163,5 @@
 			fielddisplay(obj,'damage_D','damage tensor (scalar)');
 			fielddisplay(obj,'damage_law','damage law (string) from {''initial''}');
-			if strcmpi(obj.damage_law),
+			if strcmpi(obj.damage_law,'initial'),
 				fielddisplay(obj,'damage_c1','damage parameter 1');
 				fielddisplay(obj,'damage_c2','damage parameter 2');
Index: /issm/trunk-jpl/src/m/classes/model/model.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/model/model.m	(revision 16159)
+++ /issm/trunk-jpl/src/m/classes/model/model.m	(revision 16160)
@@ -16,4 +16,5 @@
 		basalforcings    = 0;
 		materials        = 0;
+		damage           = 0;
 		friction         = 0;
 		flowequation     = 0;
@@ -187,7 +188,7 @@
 			md.materials.rheology_B=DepthAverage(md,md.materials.rheology_B);
 			md.materials.rheology_n=project2d(md,md.materials.rheology_n,1);
-			if isa(md.materials,'matdamageice')
-				md.materials.rheology_Z=DepthAverage(md,md.materials.rheology_Z);
-			end
+			
+			%damage: 
+			md.damage.D=DepthAverage(md,md.damage.D);
 
 			%special for thermal modeling:
@@ -718,7 +719,7 @@
 			md.materials.rheology_B=project3d(md,'vector',md.materials.rheology_B,'type','node');
 			md.materials.rheology_n=project3d(md,'vector',md.materials.rheology_n,'type','element');
-			if isa(md.materials,'matdamageice')
-				md.materials.rheology_Z=project3d(md,'vector',md.materials.rheology_Z,'type','node');
-			end
+				
+			%damage
+			md.damage.D=project3d(md,'vector',md.damage.D,'type','node');
 
 			%parameters
@@ -831,4 +832,5 @@
 			if isfield(structmd,'rheology_B'), md.materials.rheology_B=structmd.rheology_B; end
 			if isfield(structmd,'rheology_n'), md.materials.rheology_n=structmd.rheology_n; end
+			if isfield(structmd,'rheology_Z'), md.damage.D=(1-structmd.rheology_Z); end
 			if isfield(structmd,'spcthickness'), md.balancethickness.spcthickness=structmd.spcthickness; end
 			if isfield(structmd,'artificial_diffusivity'), md.balancethickness.stabilization=structmd.artificial_diffusivity; end
@@ -1047,4 +1049,5 @@
 			md.groundingline    = groundingline();
 			md.materials        = matice();
+			md.damage           = damage();
 			md.flowequation     = flowequation();
 			md.debug            = debug();
@@ -1079,4 +1082,5 @@
 			disp(sprintf('%19s: %-22s -- %s','basalforcings'   ,['[1x1 ' class(obj.basalforcings) ']'],'bed forcings'));
 			disp(sprintf('%19s: %-22s -- %s','materials'       ,['[1x1 ' class(obj.materials) ']'],'material properties'));
+			disp(sprintf('%19s: %-22s -- %s','damage'          ,['[1x1 ' class(obj.damage) ']'],'damage propagation laws'));
 			disp(sprintf('%19s: %-22s -- %s','friction'        ,['[1x1 ' class(obj.friction) ']'],'basal friction/drag properties'));
 			disp(sprintf('%19s: %-22s -- %s','flowequation'    ,['[1x1 ' class(obj.flowequation) ']'],'flow equations'));