Changeset 15830 for issm/trunk-jpl/src/m/mech/strainrateuncert.m

Timestamp:

08/16/13 14:50:18 (12 years ago)

Author:

cborstad

Message:

updated strain rate uncertainty calculation

File:

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

issm/trunk-jpl/src/m/mech/strainrateuncert.m

@@ -6 +6 @@
 %   The results are stored in md.results
 %
-%        'dvx' and 'dvy' are velocity errors in x and y components.  These 
-%        can either be scalars or arrays of length md.mesh.numberofvertices
+%        'dvx' and 'dvy' are velocity errors in x and y components in m/yr.  
+%   These can either be scalars or arrays of length md.mesh.numberofvertices
 %
 %   Usage:
@@ -42 +42 @@
 index=md.mesh.elements;
 summation=[1;1;1];
+dvxlist=dvx(index);
+dvylist=dvy(index);
 
 %compute nodal functions coefficients N(x,y)=alpha x + beta y +gamma
@@ -48 +50 @@
 strainrateuncert=struct('xx',[],'yy',[],'xy',[],'principalvalue1',[],'principalvalue2',[],'effectivevalue',[]);
 
-strainrateuncert.xx=dvx.*sqrt(alpha.^2*summation);
-strainrateuncert.yy=dvy.*sqrt(beta.^2*summation);
-strainrateuncert.xy=0.5*sqrt(dvx.^2.*(beta.^2*summation)+dvy.^2.*(alpha.^2*summation));
+strainrateuncert.xx=sqrt((dvxlist.*alpha).^2*summation);
+strainrateuncert.yy=sqrt((dvylist.*beta).^2*summation);
+strainrateuncert.xy=0.5*sqrt((dvxlist.*beta).^2*summation+(dvylist.*alpha).^2*summation);
 
-exx=md.results.strainrate.xx;
+exx=md.results.strainrate.xx; 
 eyy=md.results.strainrate.yy;
 exy=md.results.strainrate.xy;
-p1a=strainrateuncert.xx.*(0.5+0.25*(0.5*((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
-p2a=strainrateuncert.yy.*(0.5-0.25*(0.5*((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
-p3a=strainrateuncert.xy.*(((exx-eyy)/2).^(2)+exy.^2).^(-1./2).*exy;
-p1b=strainrateuncert.xx.*(0.5-0.25*(0.5*((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
-p2b=strainrateuncert.yy.*(0.5+0.25*(0.5*((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
-p3b=strainrateuncert.xy.*(-(((exx-eyy)/2).^(2)+exy.^2).^(-1./2).*exy);
+p1a=strainrateuncert.xx.*(0.5+0.25*(((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
+p2a=strainrateuncert.yy.*(0.5-0.25*(((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
+p3a=strainrateuncert.xy.*((((exx-eyy)/2).^2+exy.^2).^(-1./2).*exy);
+p1b=strainrateuncert.xx.*(0.5-0.25*(((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
+p2b=strainrateuncert.yy.*(0.5+0.25*(((exx-eyy)/2).^2+exy.^2).^(-1./2).*(exx-eyy));
+p3b=strainrateuncert.xy.*(-(((exx-eyy)/2).^2+exy.^2).^(-1./2).*exy);
 strainrateuncert.principalvalue1=sqrt(p1a.^2+p2a.^2+p3a.^2);
 strainrateuncert.principalvalue2=sqrt(p1b.^2+p2b.^2+p3b.^2);