Changeset 19466

Timestamp:

07/30/15 00:32:51 (10 years ago)

Author:

cborstad

Message:

CHG: working on better handling of damage-dependent mechanical properties calculation

File:

• issm/trunk-jpl/src/m/mech/mechanicalproperties.py (modified) (3 diffs)

issm/trunk-jpl/src/m/mech/mechanicalproperties.py

@@ -4 +4 @@
 from averaging import averaging
 
-def mechanicalproperties(md,vx,vy):
+def mechanicalproperties(md,vx,vy,**kwargs):
         """
         MECHANICALPROPERTIES - compute stress and strain rate for a goven velocity
@@ -29 +29 @@
         if npy.any(md.flowequation.element_equation!=2):
                 print 'Warning: the model has some non SSA elements. These will be treated like SSA elements'
+
+        #unpack kwargs
+        if 'damage' in kwargs: 
+            damage=kwargs.pop('damage')
+            if len(damage)!=md.mesh.numberofvertices:
+                raise ValueError('if damage is supplied it should be of size ' + md.mesh.numberofvertices)
+            if npy.ndim(damage)==2:
+                damage=damage.reshape(-1,)
+        else: damage=None
 
         if npy.ndim(vx)==2:
@@ -75 +84 @@
                 location=npy.nonzero(npy.logical_and(second_inv==0,power!=0))
                 nu[location]=10^18
-        elif 'matdamageice' in md.materials.__module__:
+        elif 'matdamageice' in md.materials.__module__ and damage is not None:
                 print 'computing damage-dependent properties!'
-                # FIXME might be using damage from a solution, not the initial md.damage.D
-                Zinv=npy.dot(1-md.damage.D[index-1],summation/3.).reshape(-1,)
+                Zinv=npy.dot(1-damage[index-1],summation/3.).reshape(-1,)
                 location=npy.nonzero(second_inv)
                 nu[location]=Zinv[location]*B_bar[location]/npy.power(second_inv[location],power[location])